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EXECUTIVE  SUMMARY 


Coastal  Ocean  Dynamics  Application  Radar  (CODAR)  is  a  shore-based  high-frequency  (HF) 
surface  wave  radar  system  that  can  be  used  to  measure  surface  currents  up  to  200  km  from  shore. 
Because  of  CODAR’s  ability  to  provide  real-time  surface  current  data,  the  U.  S.  Coast  Guard 
(CG)  is  investigating  the  utility  of  CODAR  in  Search  and  Rescue  (SAR)  mission  planning. 

In  2003,  the  University  of  Connecticut  (UCONN)  developed  a  Short-Term  Predictive  System 
(STPS)  that  used  CODAR  data  to  make  surface  current  predictions  24  hours  into  the  future.  The 
STPS  was  initially  developed  using  CODAR  data  from  the  standard  CODAR  system  covering 
Block  Island  Sound  (BIS).  CODAR-based  STPS  predictions  were  found  to  be  significantly  more 
accurate  than  those  developed  using  National  Oceanic  and  Atmospheric  Association  (NOAA) 
tide  data.  Based  on  this  initial  success,  additional  CODAR  development  was  pursued.  Specific 
development  initiatives  included  in  the  present  study  were  to: 

•  Extend  the  STPS  to  an  existing  long-range  CODAR  system  that  is  not  dominated  by  tidal 
flow  (i.e.,  Mid-Atlantic  Bight  (MAB).  Apply  and  evaluate  the  effect  of  geometric 
dilution  of  precision  (GDOP)  filtering  of  CODAR  data  in  the  enhanced  version  of  the 
STPS.  In  this  instance,  GDOP  filtering  was  found  to  substantially  decrease  STPS  current 
prediction  errors. 

•  Upgrade  STPS  by  adding  the  capability  to  incorporate  the  effects  of  wind  forcing  and 
evaluate  the  impact  of  including  wind  forcing.  In  addition,  the  STPS  was  enhanced  by 
the  use  of  covariance  estimates  in  the  Gauss-Markov  procedure.  Adding  forecast  winds 
improved  the  average  accuracy  of  current  forecasts  by  approximately  10  to  20  percent  in 
the  first  5  hours  of  the  forecast  in  the  MAB  region  depending  on  the  quality  of  the 
CODAR  data  used  to  make  the  prediction.  The  improvement  was  less  for  longer 
prediction  times.  This  limited  impact  of  the  wind  component  may  reflect  the  lack  of 
significant  wind  events  during  the  study  period  or  inadequacies  in  the  statistical  approach 
to  incorporate  winds. 

•  Further  investigate  the  accuracy  and  sources  of  error  in  using  CODAR  data  and  STPS  to 
determine  drift  trajectories  as  compared  to  more  conventional  methods  of  estimating 
surface  drift  (including  NOAA  tidal  data,  Mariano  climatology  data,  and  last  known 
position).  This  was  accomplished  for  both  Block  Island  Sound  and  the  Middle  Atlantic 
Bight  by  computing  separation  distance  between  model-predicted  and  actual  Self- 
Locating  Datum  Marker  Buoys  (SLDMB)  positions.  In  both  the  BIS  and  MAB,  the 
evaluation  showed  that  trajectories  calculated  from  real-time  CODAR  data  provided 
some  improvement  in  predicting  SLDMB  locations  relative  to  conventional  prediction 
approaches.  The  STPS  provided  a  measurable  improvement  over  conventional 
techniques  only  in  BIS. 

•  Investigate  methods  of  determining  and  incorporating  error  and  sub-grid  uncertainty 
terms  into  Monte  Carlo  trajectory  models  derived  from  the  comparison  of  CODAR  and 
SLDMB  trajectory  data.  In  this  part  of  the  study,  random  walk  (RW)  and  random  flight 
(RF)  methods  were  compared  for  simulating  the  transport  of  surface  drifters.  Differences 
in  RW  and  RF  search  area  size  are  small  when  comparable  input  parameters  were  used. 
Both  models  were  highly  sensitive  to  the  effect  of  uncertainties  in  dispersion. 
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•  Identify  the  current  quality  and  utility  of  HF  radar  products  for  SAR.  For  the  MAB  in 
summer,  it  was  determined  that  CODAR  improvements  did  reduce  drifting  body 
prediction  error  by  approximately  20  percent  relative  to  present  techniques. 

Improvements  in  the  processing  and  filtering  of  CODAR  returns  provided  the  largest 
benefit  in  this  study,  suggesting  that  additional  work  in  this  area  may  be  called  for.  The 
time  series  analysis  techniques  employed  in  this  study  produced  smaller  and  less 
consistent  improvements. 

•  Evaluate  the  capabilities  of  the  prototype  SAROPS  environmental  data  server  (EDS), 
which  was  extended  to  retrieve  CODAR  and  STPS  surface  current  data  and  forecasts 
from  the  MAB  CODAR  region.  This  effort  demonstrated  that  the  SAROPS  EDS  is 
capable  of  collecting  and  disseminating  CODAR  current  data  and  STPS  predictions. 

•  Conduct  an  operational  test  and  evaluation  (T&E)  at  CG  Group  Moriches,  in  which 
CODAR  and  STPS  data  were  utilized  in  actual  SAR  mission  planning.  The  T&E  showed 
that  STPS  and  CODAR  data  could  be  effectively  acquired,  interpreted,  and  used  by  CG 
SAR  mission  planners. 

CODAR  and  STPS  are  maturing  technologies  that  have  demonstrated  the  ability  to  improve  CG 
operational  planning.  Presently,  however,  CODAR  and  STPS  are  not  ready  for  operational  use 
in  SAROPS.  The  following  recommendations  are  provided  to  make  them  viable  for  operational 
use  in  SAROPS. 

•  Through  Oceans.US,  the  CG  should 

a.  Encourage  initiatives  to  improve  the  accuracy  and  coverage  of  HF  radar  arrays, 

b.  Encourage  development  of  enhancements  to  STPS  in  HF  radar  coverage  regions. 

c.  Partner  with  the  HF  radar  community  to  ensure  that  the  CG  requirements  are 
recognized. 

•  In  the  SAROPS  development  process,  G-OPR  should  include  methodologies  in  the  EDS 
to  fill  the  gaps  where  no  data  or  invalid  surface  current  data  exist. 

•  In  conjunction  with  the  Operations  System  Center,  G-OPR  should  develop  an  archive  to 
store  SLDMB  data,  and  periodically  calculate  and  store  regionally  representative 
dispersion  coefficients  from  SLDMB  trajectories.  Archived  SLDMB  data  could  also  be 
matched  with  seasonal  and  environmental  conditions  to  provide  better  historical  current 
estimates  in  areas  where  CODAR/STPS  is  not  available. 

•  G-OPR  should  support  further  research  on  the  impact  of  forecast  wind  change  on  STPS 
current  estimates. 
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1.  INTRODUCTION 


1.1.  Background 

The  U.S.  Coast  Guard  (CG)  Research  and  Development  Center  (R&DC)  conducted  research, 
development,  testing,  and  evaluation  (RDT&E)  to  evaluate  high  frequency  (HF)  surface  wave 
radar  systems  that  exploit  the  Doppler  shift  in  the  Bragg  scatter  of  the  sea  surface  to  measure 
surface  currents  in  real  time.  Several  surface  current  mapping  systems  (e.g.  Coastal  Ocean 
Dynamics  Application  Radar  (CODAR),  WEllen  Radar  (WERA),  and  Ocean  Surface  Currents 
Analyses-Real  time  (OSCAR))  have  been  developed  and  are  now  commercially  available 
(Teague,  Vesecky,  &  Fernandez,  1997).  CODAR  SeaSonde  systems  are  the  HF  surface  wave 
radar  systems  currently  used  in  the  Block  Island  Sound  (BIS)  and  Mid- Atlantic  Bight  (MAB) 
regions  covered  in  this  report.  Although  the  name  “CODAR”  is  used  throughout  this  report, 
other  brands  of  hardware  could  have  been  used  had  they  been  available  for  the  study  area. 
CODAR  systems  are  owned  and  operated  by  university  oceanographers  and  provide  surface 
current  data  for  many  areas  adjacent  to  the  U.S.  coastline  via  the  Internet.  With  recent  advances 
in  processing  speed  and  Internet-based  data  dissemination,  ocean  current  measurements  provided 
by  CODAR  may  improve  CG  mission  effectiveness  and  efficiency  in  carrying  out  search  and 
rescue  (SAR)  cases. 

In  2002,  the  R&DC  funded  a  study  (Ullman,  O'Donnell,  Edwards,  Fake,  Morschauser,  Sprague, 
Allen,  &  Krenzien,  2003)  with  participation  by  the  University  of  Connecticut  (UCONN)  and 
University  of  Rhode  Island  (URI)  to  investigate  more  closely  the  application  of  data  from  the 
Block  Island  Sound  (BIS)  CODAR  system  in  Search  and  Rescue  (SAR)  planning.  The  BIS 
CODAR  system  is  a  standard-range  25 -MHz  CODAR  system,  owned  and  operated  by  URI  and 
UCONN.  It  operates  from  three  stations:  Montauk  Point,  NY;  Southeast  Light  House,  Block 
Island,  RI;  and  Misquamicut,  RI.  The  locations  of  the  BIS  CODAR  transmitters  and  the 
approximate  operating  ranges  of  each  transmitter  are  shown  in  Figure  1-1.  Using  direction-  and 
range-finding  algorithms,  each  station  measures  radial  surface  current  velocities  every 
10  minutes.  These  radial  velocities  are  filtered  and  used  to  calculate  the  surface  current 
velocities  for  the  BIS,  a  region  approximately  40  x  40  kilometers  (km),  with  grid  cells  every 
1.5  km.  The  depth  of  the  water  column  over  which  currents  are  measured  is  dependent  on  the 
transmitter  frequency:  at  25-MHz,  surface  currents  in  the  upper  0.5  meters  of  the  water  column 
are  measured.  Theoretically,  the  current  field  is  resolved  in  areas  where  the  arcs  of  at  least  two 
of  the  transmitters  intersect.  Coverage  is  lost,  however,  in  the  shadow  of  a  land  mass,  such  as 
Block  Island,  or  in  the  line  of  sight  between  two  stations,  or  at  the  outer  range  of  each 
transmitter.  The  system  produces  surface  current  maps  at  one  hour  intervals.  A  map  depicting 
surface  current  typical  coverage  of  the  BIS  system  is  shown  in  Figure  1-2. 

Ullman,  et  al,  (2003)  focused  on  the  following  three  tasks.  (1)  develop  a  Short-term  Predictive 
System  (STPS)  based  upon  surface  current  fields  estimated  from  CODAR  data  in  the  Block 
Island  Sound  (BIS)  region;  (2)  develop  an  environmental  data  server  (EDS)  Web  site  for 
collecting  and  disseminating  CODAR  data  fields  and  forecast  fields  from  the  STPS;  and 
(3)  compare  the  actual  trajectories  of  surface  drifters  with  predicted  trajectories  based  upon 
CODAR,  STPS,  and  the  nearest  National  Oceanographic  and  Atmospheric  Administration 
(NO  A  A)  tidal  current  prediction  stations. 
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Figure  1-1.  The  coastline  and  bathymetry  of  Block  Island  Sound  with  the  location  and 
approximate  observation  areas  of  three  CODAR  sites. 
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Figure  1-2.  The  typical  spatial  coverage  of  the  BIS  CODAR  system  operated  by  URI  and 
UCONN. 


1-2 


In  the  present  project,  the  features  of  the  BIS  CODAR  system  is  extended  to  the  CODAR 
network  for  the  Middle  Atlantic  Bight  (MAB).  The  MAB  CODAR  system  is  a  long-range 
system  consisting  of  ten  stations  operated  by  Rutgers  University.  Six  of  these  stations  are  lower- 
resolution,  long-range  sites.  With  an  approximate  5-MHz  operating  frequency,  these  sites 
measure  surface  currents  within  the  upper  2.4  meters  of  the  water  column  (Stewart  and 
Joy,  1974).  Typical  spatial  resolutions  are  on  the  order  of  6  km  (grid  cell  size)  with  maximum 
ranges  exceeding  200  km.  A  cluster  of  four  sites  along  the  coast  of  New  Jersey  extending 
between  Wildwood  and  Sandy  Hook  form  the  MAB  CODAR  network  for  this  study,  which 
provides  hourly  surface  current  maps  over  an  area  approximately  150  x  100  km  off  the  New 
Jersey  coast.  These  sites  use  global  positioning  system  (GPS)  time  synchronization  so  that  each 
site  operates  at  the  same  frequency  and  is  bi-statically  linked  to  the  other  sites  in  the  network. 
The  CODAR  system  samples  data  continuously  and  outputs  a  4-hour  running  average  every 
hour.  Processing  of  the  CODAR  data  performed  by  Rutgers  University  is  further  outlined  in 
Appendix  A.  Figure  1-3  shows  the  MAB  CODAR  coverage  for  July/ August  2004.  The  black 
line  in  the  figure  bounds  the  50  percent  contour  coverage  for  the  length  of  the  deployment,  i.e., 
grid  points  inside  the  black  line  had  returns  more  than  50  percent  of  the  time;  grid  points  outside 
of  the  black  line  had  returns  less  than  50  percent  of  the  time.  The  50  percent  contour  coverage  is 
based  on  data  from  the  27  July  to  3 1  August  2004  deployment  period. 


Longitude  (Degrees:  Minutes) 

Figure  1-3.  MAB  CODAR  coverage  for  July/ August  2004. 


In  a  separate  effort,  the  CG  has  been  developing  a  system  for  search  planners,  called  Search  and 
Rescue  Optimal  Planning  System  (SAROPS).  SAROPS  is  intended  to  replace  both  of  the  CG’s 
presently  fielded  search  planning  tools,  Command  and  Control  Personal  Computer  (C2PC)  SAR 
Tools  and  Computer-Assisted  Search  Planning  (CASP).  A  prototype  graphical  user  interface 
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(GUI)  has  been  provided  by  the  SAROPS  developer  as  part  of  this  study  to  see  how  SAR 
controllers  could  utilize  the  current  fields  provided  by  CODAR  in  SAROPS  to  plan  search 
efforts.  Procedures  and  results  are  presented  in  Section  6  of  this  report. 

1.2.  Project  description 

The  goal  of  this  project  is  to  resolve  the  potential  benefits  to  SAR  Planning  of  adapting  the  STPS 
to  a  larger  and  more  generalized  CODAR  network.  The  elements  of  this  project  are  (1)  modify 
the  STPS  developed  for  the  tidally-dominated  BIS  to  the  larger  MAB  where  circulation  is 
dominated  by  winds  and  other  processes;  (2)  develop  and  exercise  new  analysis  procedures  to 
incorporate  wind  forcing  into  STPS;  (3)  identify  the  processes  or  factors  that  affect  the 
performance  of  the  CODAR-based  surface  current  predictions  in  the  MAB;  and  finally  (4) 
conduct  an  operational  test  of  the  system  to  identify  its  skill  and  to  identify  the  direction  and 
potential  benefits  of  future  work  in  this  area. 
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2.  MODIFICATIONS  TO  STPS  FOR  THE  MAB  APPLICATION 


The  STPS  developed  in  Ullman  et  al,  (2003)  is  based  on  the  decomposition  of  the  CODAR 
measured  currents  into  tidal  and  non-tidal  components  (i.e.,  wind-driven  motion  and  motion 
associated  with  ocean  density  variations  and  very  low  frequency  waves).  At  each  grid  point, 
tidal  currents  are  predicted  using  tidal  constituents  derived  from  a  harmonic  analysis  of  a  one- 
month  record  of  CODAR  surface  currents.  A  Gauss-Markov  estimator  is  used  to  predict  the  non- 
tidal  component  of  currents.  The  current  prediction  is  the  sum  of  the  predicted  tidal  and  non- 
tidal  currents. 

The  Gauss-Markov  estimator,  based  on  auto-  and  cross-covariance  functions,  is  used  to  predict 
the  non-tidal  component  of  current  (Ullman,  O'Donnell,  Kohut,  Fake,  Allen,  2005).  The 
covariance  function  is  a  series  that  measures  how  similar  a  time  series  signal  is  to  a  time-shifted 
version  of  itself  or  another  series  of  values.  For  example,  in  a  tidally-dominated  area,  the 
autocovariance  series  will  have  a  high  value  at  time  lags  corresponding  to  tidal  periodicities  (e.g. 
12.4  hours).  The  Gauss-Markov  estimator  calculates  filter  weights  that  also  are  a  function  of 
time  lag,  which  are  then  used  to  generate  a  prediction  of  future  values  based  on  a  memory  effect 
created  from  previous  values  through  the  covariance  function.  Predictions  are  produced  on  an 
hourly  basis  and  extend  for  25  hours  into  the  future. 

The  Gauss-Markov  estimator  is  used  in  two  ways  in  this  study.  In  the  first  application,  it  derives 
a  series  of  weights  that  are  calculated  from  the  autocovariance  of  previous  current  data  that  are 
then  applied  to  recent  values  to  predict  future  currents  (i.e.  for  the  next  25  hours).  In  the  second 
application,  the  cross-covariance  between  previous  surface  wind  measurements  and  surface 
currents  derived  from  CODAR  will  be  used  to  estimate  future  currents. 

The  STPS  developed  for  BIS  did  not  explicitly  account  for  wind  forcing  in  the  BIS  because  its 
circulation  is  considered  to  be  tidally  dominated.  In  exercising  the  STPS  and  comparing  the 
results  to  SLDMB  trajectories  in  the  Middle  Atlantic  Bight,  it  seemed  likely  that  a  second 
generation  model  that  included  wind  forcing  could  improve  forecast  accuracy  in  areas  where 
there  is  a  significant  wind-driven  surface  current  component. 

The  STPS  algorithm  was  extended  to  accommodate  data  from  the  MAB  between  southern  Long 
Island  and  the  mouth  of  Delaware  Bay.  The  following  modifications  were  implemented  and 
tested  in  STPS: 

•  Covariance  calculated  for  each  grid  point,  instead  of  a  region-wide  average,  was  used  in 
the  Gauss-Markov  procedure  to  account  for  local  geographic  and  bathymetric  effects. 

•  Covariance  was  calculated  between  the  east  and  north  velocity  components  of  current  at 
each  grid  point  to  account  for  the  Coriolis  effect. 

•  Improvements  were  incorporated  into  the  MAB  CODAR  processing  algorithm,  including 
geometric  dilution  of  precision  (GDOP)  and  percentage  return  filtering. 

•  Measured  wind  data  from  a  moored  buoy  in  the  vicinity  of  the  CODAR  coverage  area 
was  used  to  calculate  the  covariance  function  with  the  CODAR  data.  The  wind  data  were 
used  in  two  ways.  First,  wind  observations  prior  to  each  time  step  were  used  to  predict 
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future  surface  currents  and  second,  “predicted”  winds,  (actually  wind  observations)  were 
used  to  predict  future  surface  currents. 

The  incorporation  of  local  covariance  and  east  and  north  velocity  component  covariance  into  the 
Gauss-Markov  estimator  (items  1  and  2  above)  did  not  provide  significant  improvements  to  the 
model’s  predictive  skill.  These  improvements  were  retained  in  the  STPS,  but  are  not  discussed 
further.  Details  of  the  STPS  enhancements  in  the  BIS  region  are  described  in  detail  in 
Appendix  B. 

2.1.  Enhancements  adopted  from  the  MAB  algorithm 

To  further  refine  the  CODAR  data  for  the  present  project,  Rutgers  University  implemented  a 
new  algorithm  for  combining  CODAR  radials  that  incorporated  two  improvements:  a  50  percent 
return  mask  and  GDOP  filtering  (Kohut,  Roarty,  Glenn,  2004).  Each  hour,  the  available  radial 
vector  maps  were  combined  into  a  single  total  vector  map.  All  radial  component  vectors  within 
10  km  were  used  to  calculate  the  vector  at  each  grid  point,  but  a  total  vector  was  generated  only 
if  at  least  three  radial  vectors  from  at  least  two  remote  sites  were  used  in  the  combination.  The 
uncertainty  of  the  combined  totals  was  separated  into  radial  vector  uncertainty  based  on  the 
number  of  radials  in  the  4-hour  average,  signal  to  noise  ratio  (SNR),  and  geometric  dilution  of 
precision  (GDOP)  (e.g.  see  Riddles,  2003).  Scalar  values  (Chapman  and  Graber,  1997)  were 
calculated  based  on  the  intersection  angles  of  radial  vectors  within  10  km  of  each  grid  point  and 
compared  with  a  scalar  threshold.  The  expected  error  in  the  surface  current  vector  due  to  GDOP 
calculated  from  two  radials  increases  as  their  intersection  angle  diverges  from  the  orthogonal. 

For  this  project,  the  scalar  GDOP  threshold  was  calculated  to  be  1.25,  corresponding  to 
intersection  angles  less  than  30  degrees  or  greater  than  150  degrees.  These  filters  reduced  the 
size  of  the  CODAR  coverage  area,  but  they  improved  the  quality  of  the  data. 

Figure  2-1  shows  MAB  CODAR  coverage  when  the  filtering  algorithm  was  used  for 
July/ August  2004.  The  black  line  in  the  figure  bounds  the  50  percent  contour  coverage  for  the 
length  of  the  deployment.  The  50  percent  contour  coverage  is  based  on  data  from  the  27  July  to 
3 1  August  2004  deployment  period.  When  the  GDOP  and  percent  return  constraints  were  used 
to  filter  the  set  of  CODAR  data  used  to  generate  current  predictions  (referred  to  below  as  the 
“New  Gauss-Markov”  predictor),  the  root-mean- squared  (RMS)  error  of  the  predictions  was 
reduced  by  approximately  20  percent. 
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Figure  2-1.  MAB  CODAR  coverage  for  July /August  2004,  showing  percent  radar  return 
percentages  (color  coded)  and  using  the  GDOP-filtered  data.  (Note:  The 
approximate  location  of  NOAA  buoy  44025  is  at  the  top  of  the  50  percent 
coverage  area.) 


2.2.  Incorporation  of  NOAA  buoy  data  into  STPS 

The  first  attempt  to  incorporate  wind  data  into  STPS  involved  the  use  of  previous  wind 
observations  to  construct  predictions  of  near  future  non-tidal  residual  currents.  Wind  data  from 
the  closest  operating  NOAA  buoy  (#44025)  were  used  for  this  exercise.  We  chose  a  seven  day 
period  during  August  2004  when  winds  were  highest  (~5  m/s)  to  maximize  the  influence  of 
wind.  First,  previous  winds  and  CODAR  current  fields  were  used  to  establish  filter  weights  for 
the  Gauss-Markov  estimator.  Next,  at  each  time  step,  the  wind-driven  residual  current  was 
predicted  for  the  upcoming  24  hour  period.  This  prediction  was  made  with  the  Gauss-Markov 
estimator  using  observed  winds  and  CODAR  current  fields  from  the  previous  24  hour  period. 
For  the  purpose  of  identifying  the  degree  of  improvement,  an  initial  run  was  made  where  winds 
were  set  to  zero,  such  that  only  the  basic  features  of  the  STPS,  tides  plus  a  residual  current 
estimated  using  the  Gauss-Markov  estimator,  were  in  play. 

The  improvement  provided  by  using  recent  buoy  winds  with  the  basic  STPS  for  August  2004  is 
summarized  in  Figure  2-2  by  plotting  the  spatially-averaged  RMS  differences  between  the  STPS 
and  CODAR  surface  current  components  versus  forecast  time  for  the  August  study  period.  The 
“New  Gauss-Markov”  titles  in  the  figure  signify  that  the  CODAR  return  processing  algorithms 
described  in  Section  2-1  were  being  used.  The  no-wind  case  differences  are  shown  as  a  function 
of  forecast  time  for  the  east  (u-component,  blue  solid  line)  and  north  (v-component,  red  solid 
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line)  current  components.  Where  the  wind-driven  residual  is  not  incorporated,  the  RMS 
differences  start  at  9  cm/s  at  the  t=0  forecast  and  increase,  then  level  off  at  12  cm/s  for  the  east 
component  and  13  cm/s  for  the  north  component  at  approximately  5  hours  forecast  time.  The 
dashed  blue  line  (u-component)  and  dashed  red  line  (v-component)  show  the  corresponding 
RMS  difference  plots  incorporating  the  wind-driven  residual  current  based  on  observed  wind. 

The  reductions  in  RMS  error  vs.  time  plots  at  the  bottom  of  the  figure  show  that  a  reduction  1  to 
2  cm/s  (20  percent)  was  realized  for  the  present  time  analysis,  but  the  reduction  dropped  off 
rapidly  to  lcm/s  or  less  at  forecast  times  of  a  few  hours  or  more. 


Time  Since  Prediction  (h) 

Figure  2-2.  The  RMS  difference  (cm/s)  between  area-averaged  MAB  CODAR  east  (u-blue 
lines)  and  north  (v-red  lines)  components  of  surface  current  for  August  2004  as  a 
function  of  current  forecast  period  (hrs)  using  measured  wind  (dashed  lines)  and 
no  wind  forcing  (solid  lines).  The  reductions  (Au,  Av)  in  RMS  error  at  the  bottom 
of  the  figure  represent  the  error  reduction  realized  by  the  use  of  buoy  winds. 

2.3.  Use  of  “forecast”  winds  derived  from  NOAA  buoy  data 

An  experiment  was  performed  in  which  observed  winds  from  NDBC  buoy  #44025  were 
included  as  “forecast”  winds  during  the  same  period  in  August  2004.  In  this  test  winds  from 
times  following  each  present  time  (i.e.  as  simulated  predictions  of  future  winds)  were  included  in 
the  Gauss-Markov  estimator.  In  making  predictions  for  planning  an  actual  SAR  case,  these 
“future  wind  values”  would  of  course  be  replaced  by  model  predictions.  In  a  sense,  this  provided 
a  “perfect  prediction”  of  future  winds  for  the  time  period  considered  to  establish  an  upper  limit 
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of  the  potential  value  of  including  wind  forecasts.  At  each  time  step,  the  “New  Gauss-Markov 
Forecast”  approach  was  used  to  construct  an  estimate  of  the  surface  current  based  on  measured 
buoy  winds  during  the  48  hour  period  centered  on  the  present  time  (e.g.  t  =  0  in  Figure  2-3).  As 
before,  the  basic  features  of  the  STPS  (Gauss-Markov  tidal  and  residual  currents)  were  included 
in  the  prediction.  The  improvement  provided  by  using  present  buoy  winds  relative  to  the  basic 
STPS  using  only  observations  (made  prior  to  t=0)  is  shown  by  plotting  RMS  difference  as  a 
function  of  forecast  time  in  Figure  2-3  for  the  August  2004  time  period.  The  area-averaged  east 
(blue)  and  north  (red)  RMS  velocity  component  errors  are  plotted  as  a  function  of  the  forecast 
period  (hrs).  The  solid  lines  represent  predictions  based  on  measured  and  forecast  winds;  dashed 
lines  correspond  to  predictions  based  only  on  measured  winds.  The  results  suggest  that  the  use 
of  wind  forecasts  produces  a  negligible  improvement  in  prediction  skill.  The  low  magnitude  and 
variability  of  the  wind  during  the  test  period  may  not  have  produced  a  significant  response  from 
the  model  because  surface  currents  were  not  well  correlated  with  local  winds  during  that  time. 
Other  explanations  may  be  that  surface  currents  are  forced  by  other  processes,  or  that  another 
analysis  approach  to  quantifying  the  wind-driven  circulation  is  called  for. 
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Time  Since  Prediction  (h) 

Figure  2-3.  The  evolution  of  the  area-averaged  RMS  difference  between  the  forecast  current 

components  using  the  enhanced  STPS  algorithm  with  measured  wind  (dashed  lines) 
and  with  measured  and  forecast  wind  (solid  lines)  in  the  MAB  area.  The  reductions 
(Au,  Av)  in  RMS  error  at  the  bottom  of  the  figure  represent  the  error  reduction 
realized  by  the  use  of  buoy  winds. 
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2.4.  Summary  of  STPS  extension  efforts 


A  series  of  changes  were  made  to  the  STPS  as  it  was  adapted  from  the  BIS  area  to  the  MAB 
area.  The  changes  fell  into  two  categories:  modifications  to  the  processing  of  the  CODAR  data 
to  enhance  data  quality  prior  to  use  by  STPS,  and  application  of  new  techniques  to  predict 
surface  currents  incorporating  a  wind-driven  component  based  on  actual  observed  and  forecast 
winds.  The  exercises  determined  that  changes  to  the  processing  of  CODAR  returns  had  the 
largest  influence  on  the  agreement  of  subsequent  CODAR-derived  surface  current  data  with 
STPS  predictions.  The  use  of  currents  derived  from  the  GDOP  and  percent  return  filters  reduced 
RMS  errors  by  approximately  20  percent. 

The  introduction  of  wind  data  into  the  Gauss-Markov  estimator  as  incorporated  into  the  MAB 
STPS  algorithm  provided  a  small  improvement  in  the  RMS  differences  between  predictions  and 
future  CODAR  data.  The  reductions  in  RMS  prediction  error  varied  between  2  cm/s  at  the  one 
hour  forecast  time  to  less  than  1  cm/s  at  prediction  times  of  a  few  hours,  or  approximately  10 
percent. 

Potential  reasons  for  the  minimal  contribution  of  the  wind-driven  component  may  be  that  the  low 
magnitude  and  variability  of  the  wind  during  the  test  period  may  not  have  produced  a  significant 
response  from  the  model  because  surface  currents  were  not  well  correlated  with  local  winds 
during  that  time.  Other  reasons  may  be  that  the  statistical  approach  chosen  did  not  represent  the 
effect  of  wind.  Additional  testing  and  a  different  conceptual  approach  are  needed,  particularly 
for  the  0  to  6  hour  interval,  a  time  period  of  significant  operational  interest  to  the  CG. 

In  summary,  the  greatest  improvement  in  the  STPS  performance  was  due  to  improved  CODAR 
data  processing  than  to  adding  the  effect  of  forecast  winds  to  the  STPS. 
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3.  EVALUATION  OF  PREDICTIVE  METHODS  IN  ESTIMATING 
DRIFTER  TRAJECTORIES 

Through  this  study  and  previous  work,  investigators  had  access  to  actual  SLDMB  trajectory 
information  in  the  Mid-Atlantic  Bight  (MAB)  and  Block  Island  Sound  (BIS)  dating  from 
December  2002  through  September  2004.  These  trajectories  were  used  to  evaluate  the 
performance  of  the  extended  STPS  with  wind  forcing,  and  their  underlying  CODAR  surface- 
current  maps  with  improved  processing,  and  compare  the  CODAR/STPS  performance  against 
traditional  CG  search  planning  methods  for  determining  surface  currents. 

3.1.  Mid- Atlantic  bight  SLDMB  deployments 

Figure  3-1  and  Figure  3-2  show  the  trajectories  of  SLDMBs  launched  in  the  MAB  study  area. 
The  thin  dark  line  shows  the  area  for  which  CODAR  coverage  was  available  more  than 
90  percent  of  the  time.  The  thick  lines  represent  the  drifter  trajectories. 

Drifter  Posilions  Between  27-Jyl-20D4  1250:00  And  D1 -Sep-2004  14:00:00 


Figure  3-1.  Drifter  trajectories  in  the  MAB  study  area  in  July-Sept  2004  (Black  dots  represent 
the  start  point  for  each  trajectory). 
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Drifter  Positions  Between  27 -Mar-2 003  03:00:00  And  1  &-Apr-2003  17:00:00 


Longitude  (°W) 


Figure  3-2.  Trajectories  in  the  MAB  study  area  in  March-April  2003  (Black  dots  represent  the 
start  point  for  each  trajectory). 


The  drifter  trajectories  were  divided  into  24-hour  segments,  creating  approximately 
200  independent  sample  tracks.  The  frequency  distribution  of  separation  distances  (position 
errors)  between  drifter  location  and  positions  computed  using  the  CODAR  current  estimates  at  3, 
6,  12,  and  24  hours  is  shown  in  the  left  column  of  Figure  3-3.  The  second  column  shows  the 
same  information  with  the  drifter  locations  predicted  using  the  STPS  currents. 

Measured  wind  data  from  NDBC  buoy  #44025  were  used  in  the  STPS  predictions  to  incorporate 
a  wind-driven  component  as  described  in  Section  2.2.  Because  the  CG  currently  uses  a 
climatological  mean  current  specified  at  1 -degree  intervals  to  estimate  the  drift  of  search  targets 
in  most  of  this  area,  the  third  column  shows  the  separation  distances  between  the  drifter  location 
and  the  positions  predicted  using  the  Mariano  climatology  (Mariano,  Ryan,  Perkins,  Smithers, 
1995).  Finally,  the  LKP  column  shows  the  distribution  of  separation  distances  between  the 
actual  position  of  the  drifter  and  its  initial  location  (i.e.,  if  we  assumed  that  the  object  did  not 
drift  from  its  initial  position  over  the  24-hour  period).  At  3  and  6  hours,  the  distributions  of  the 
separations  in  the  left  two  columns  (CODAR  and  STPS)  are  slightly  more  skewed  to  lower 
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values  than  the  right  two  columns,  indicating  that  these  two  approaches  are  somewhat  more 
accurate. 
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Figure  3-3.  Histograms  of  the  separation  between  the  predicted  and  actual  locations  of  SLDMBs 
in  the  MAB  study  area,  shown  as  the  number  of  endpoints  (N)  vs  separation 
distance  (km)  at  3,  6,  12,  and  24  hours.  From  left  to  right,  the  predicted  positions 
are  obtained  by  integrating  observed  CODAR  velocities,  STPS-predicted  currents, 
velocities  estimated  using  the  Mariano  climatology  data,  and  the  LKP. 
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Distances  corresponding  to  the  95th  percentile  of  the  separation  distance  distributions  at  each 
time  step  are  shown  in  Figure  3-4.  The  differences  are  relatively  similar  in  this  case.  At  24 
hours,  the  difference  between  the  four  methods  is  less  than  5  km.  As  shown  in  Figure  3-4, 
CODAR  and  STPS  have  slightly  smaller  separation  distances  than  the  Mariano  climatology  or 
LKP  estimates,  but  the  difference  is  not  significant. 


Figure  3-4.  The  evolution  of  the  95th  percentile  value  of  separation  distance  between  drifters 

and  simulated  trajectories  in  the  MAB  based  on  CODAR  data  (green  line),  the  STPS 
algorithm  (blue  line),  LKP  (black  line),  and  Mariano  climatology  (red  line). 


To  provide  a  more  realistic  comparison  of  the  CODAR  and  STPS  predictions  with  the  more 
conventional  predictions,  data  from  the  MAB  area  were  reprocessed  to  exclude  data  outside  of 
the  90  percent  availability  line  (i.e.,  excluding  lower  percent  return  CODAR  areas).  The  results 
are  provided  in  Figure  3-5  and  Figure  3-6.  The  data  are  similar  to  those  without  screening  except 
for  CODAR.  The  CODAR  errors  decrease  approximately  20  to  25  percent  when  the  lower- 
quality  data  are  excluded,  and  STPS  errors  decrease  by  roughly  10  percent.  The  outcome 
illustrates  the  value  gained  by  preprocessing  and  use  of  higher  quality  CODAR  returns 
implemented  in  the  MAB  processing  algorithm,  in  improving  the  performance  of  these 
approaches  for  predicting  drift  body  motion. 


3-4 


Mariano 


100 


50 


50 

60 

40 

20 


-  A  t  -  1 2h 

At  -  I2h 

At  -  12h 

At  -  I2h 

L 

L 

L 

aJ 

50 

40 

30 

20 

1D| 


A  t  -  24h 

A  t  -  24h 

A  t  -  24h; 

1 - ~~\ 

A  t  -  24h 

L 

L 

JL. 

4. 

3  10  20  30  40 
separation  (km) 


3  10  20  30  40 
separation  (km) 


J  10  20  30  40 
separation  (km) 


0  10  20  30  40 
separation  (km) 


Figure  3-5.  Histograms  of  the  separation  between  the  predicted  and  actual  locations  of  SLDMBs 
in  the  MAB  study  area,  shown  as  the  number  of  endpoints  (N)  vs  separation  distance 
(km)  at  3,  6,  12,  and  24  hours.  From  left  to  right,  the  predicted  positions  are 
obtained  by  integrating  observed  CODAR  velocities,  STPS-predicted  currents, 
velocities  estimated  using  the  Mariano  climatology,  and  the  LKP.  These  results 
were  obtained  using  only  drifter  trajectories  that  were  always  inside  the  CODAR  90 
percent  availability  line. 
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Figure  3-6.  The  evolution  of  the  95th  percentile  separation  distance  between  drifters  and 

simulated  trajectories  in  the  MAB  based  on  CODAR  data  (green  line),  the  STPS 
predictions  (blue  line),  LKP  (black  line),  and  Mariano  climatology  (red  line). 
These  results  were  obtained  using  only  drifter  trajectories  that  were  always 
inside  the  CODAR  90  percent  availability  line. 

3.2.  Block  Island  Sound  SLDMB  deployments 

The  four  predictive  approaches  were  applied  to  the  BIS  for  the  purpose  of  comparing  their 
performance  in  a  tidally-dominated  region.  The  SLDMB  deployments  also  included  trajectories 
inside  the  area  covered  by  the  BIS  CODAR  system.  Figure  3-7  shows  the  locations  of  the 
release  points  (red  dots)  of  the  drifters  and  the  subsequent  tracks  (blue).  Two  obvious  features  of 
these  trajectories  deserve  comment. 

•  The  drifters  launched  in  the  northern  half  of  the  BIS  moved  to  the  east  and  exit  the  sound 
between  Block  Island  and  the  Rhode  Island  shore.  In  contrast,  those  released  in  the 
southern  half  moved  almost  directly  to  the  south  between  Montauk  and  Block  Island  into 
the  Atlantic.  Thus,  knowing  the  initial  position  and  the  spatial  structure  of  the  circulation 
is  critical  to  knowing  where  drifters  will  go. 

•  Within  BIS,  large  reversals  occur  in  the  trajectories  where  the  tidal  currents  are  strong. 

The  drifter  trajectories  were  divided  into  24-hour  segments  overlapped  by  12  hours,  which 
resulted  in  a  total  of  140  independent  sample  tracks  within  the  area  monitored  by  the  BIS 
CODAR  system.  The  modified  STPS  was  used  to  make  24-hour  trajectory  predictions  to 
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compare  with  the  subsequent  SLDMB  tracks  along  with  predictions  based  on  actual  CODAR 
data  (without  the  MAB  pre-processing  steps)  and  the  NOAA  tidal  current  forecasts.  Position 
difference  from  the  LKP  was  used  as  a  reference.  In  each  case,  the  observed  (SLDMB)  and 
predicted  (CODAR,  STPS,  NOAA  Tide,  and  LKP)  velocity  fields  were  converted  to  trajectories 
using  simple  predictor-corrector  integration  combined  with  bi-linear  interpolation  of  point 
velocity  data  at  nearby  locations  to  the  locations  of  the  drifters. 


-7J  -71.9  -719  -71,7  -71.6  -71.5  -714  -71.3 

Lro.fW) 


Figure  3-7.  Deployment  locations  of  SLDMBs  (red  dots)  and  their  subsequent  trajectories  (blue 
lines)  in  BIS. 


The  distributions  of  position  errors  as  represented  by  separation  distances  between  observed  and 
predicted  drifter  locations  are  shown  at  3,  6,  12,  and  24  hour  forecast  times  in  Figure  3-8. 
Separation  distances  between  actual  SLDMB  locations  and  CODAR  predicted  locations  are 
shown  in  the  left  column.  The  second  column  shows  the  separation  distances  between  actual 
SLDMB  locations  and  STPS  predicted  locations.  Because  the  CG  currently  uses  NOAA  tidal 
current  station  information  to  estimate  the  drift  of  search  targets  in  most  of  this  area,  the  third 
column  shows  the  separation  distances  between  the  SLDMB  locations  and  predictions  made 
using  the  NOAA  information.  (Note:  data  from  the  NOAA  tidal  station  nearest  to  the  position  at 
the  beginning  time  segment  was  applied  for  the  duration  of  the  24-hour  trajectory.)  Finally,  to 
evaluate  the  assumption  that  an  object  did  not  drift  from  its  last  known  position  (LKP),  the  right 
column  shows  the  separation  distances  between  the  actual  position  of  the  drifter  and  its  initial 
location  estimate  (LKP). 
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Figure  3-8.  Histograms  of  the  separation  between  the  predicted  and  actual  locations  of  SLDMBs 
in  the  MAB  study  area,  shown  as  the  number  of  endpoints  (N)  vs.  separation 
distance  (km)  at  3,  6,  12,  and  24  hours.  From  left  to  right,  the  predicted  positions 
are  obtained  by  integrating  observed  CODAR  velocities,  STPS-predicted  currents, 
velocities  estimated  using  the  NOAA  tides,  and  the  LKP. 


The  NOAA  tide  and  LKP  separation  histograms  are  more  widely  distributed,  particularly  as  the 
prediction  time  increases,  indicating  that  position  errors  will  be  larger  for  these  approaches.  The 
CODAR  and  STPS  separation  distance  distributions  show  smaller  separation  distances  (position 
errors),  particularly  within  the  first  6  hours,  indicating  that  these  approaches  have  higher 
predictive  skill.  At  3  and  6  hours,  the  separation  distance  is  0  to  2  km  with  almost  all  values  less 
than  5  km  for  both  the  CODAR  and  STPS  distributions.  At  12  and  24  hours,  the  modal  values 
increase  to  4  km  and  the  distribution  of  the  STPS  separations  is  obviously  wider  than  that  in  the 
first  column,  reflecting  the  decrease  in  the  predictive  skill  with  increasing  forecast  time  for 
STPS. 
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Figure  3-9  provides  another  view  of  the  separation  by  showing  the  increase  in  the  95th  percentile 
of  each  distribution  in  Figure  3-8  as  a  function  of  prediction  time.  The  CODAR  analysis 
trajectories  and  STPS  predictions  produce  nearly  equivalent  separation  distance  (position  errors) 
between  0  and  5  hours  with  an  approximate  RMS  separation  error  of  5  km  at  5  hours.  On  the 
other  hand,  the  separation  based  on  the  NOAA  tidal  current  predictions  was  approximately 
20  km  at  5  hours,  a  factor  of  4  higher.  Of  the  four  methods,  the  NOAA  tidal  current  predictions 
provided  the  least  accurate  estimate  of  drifter  location.  This  result  is  partially  due  to  the  distance 
and  position  of  the  drifters  relative  to  the  nearest  NOAA  tidal  station  and  because  the  NOAA 
tidal  stations  used  for  estimates  were  not  updated  during  the  drift.  The  no-motion  drifter 
displacement  (LKP)  error  is  lower  than  the  NOAA  tides  prediction  with  a  separation  of 
approximately  7.5  km  at  5  hours  but  is  less  accurate  than  the  STPS  forecast.  Beyond  the  5-hour 
prediction  time,  the  STPS  prediction  diverges  from  the  CODAR-derived  positions,  and  is  much 
larger  at  the  24-hour  prediction  time.  The  separation  distance  in  SLDMB  positions  associated 
with  the  use  of  CODAR  data  is  approximately  16  km  after  a  24-hour  period. 

Bl  Region:  95th  Percentile  Separation 


Figure  3-9.  The  evolution  of  the  95th  percentile  separation  distance  between  drifters  and 
simulated  trajectories  in  the  BIS  based  on  CODAR  data,  the  STPS  algorithm, 
last  known  position  (LKP),  and  NOAA  tidal  current  predictions. 

3.3.  Summary  of  results 

In  both  BIS  and  MAB,  the  CODAR  systems  provided  the  most  accurate  estimates  of  the  four 
methods  selected.  This  result  validates  the  usefulness  of  CODAR,  but  is  tempered  somewhat  by 
the  fact  that  CODAR  is  a  real-time  observation  and  cannot  be  used  as  a  predictive  tool  for  these 
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exercises.  In  both  areas,  the  STPS  produced  larger  separation  distances  than  the  CODAR 
analysis  velocities  and  provided  only  a  modest  improvement  in  prediction  accuracy  over 
conventional  methods.  Use  of  the  Mariano  currents  (standard  CG  procedure)  in  the  MAB 
provided  a  prediction  of  SAR  object  drift  that  was  functionally  identical  to  the  STPS.  This  result 
is  reasonable  because  winds  were  relatively  light  over  the  period  of  the  test.  For  the  MAB, 
trajectory  errors  derived  from  CODAR  data  improved  the  skill  of  the  prediction,  reducing  the 
separation  error  by  20  to  25  percent  when  CODAR  data  outside  of  the  90  percent  availability  line 
were  not  used  for  trajectory  calculations.  In  BIS,  trajectory  predictions  based  on  the  NOAA  tidal 
current  predictions  produced  the  largest  prediction  errors  of  any  method.  This  is  caused  by  the 
spatial  variability  of  tidal  currents  in  that  area  and  because  current  estimates  were  used  only  from 
the  closest  operational  NOAA  station  to  the  study  area  at  the  start  of  the  period  of  drift.  The 
signal  processing  and  use  of  wind  data  in  STPS  did  not  improve  or  change  its  predictive  skill  in 
BIS. 

CODAR  therefore  represented  a  modest  improvement  over  the  present  Mariano  current  data  tool 
in  this  area.  The  use  of  STPS  produced  only  a  minor  improvement  over  those  derived  from  the 
LKP  and  NOAA  tidal  current  predictions.  The  incorporation  of  the  MAB  CODAR  filtering 
algorithms  produced  the  largest  improvement  in  system  accuracy. 
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4.  OTHER  FACTORS  AFFECTING  USABILITY  OF  CODAR  AND  STPS 
FOR  SEARCH  PLANNING 

This  section  focused  on  the  characterization  of  uncertainties  in  the  surface  current  field  using 
CODAR  and  STPS,  and  the  subsequent  incorporation  of  these  parameters  into  trajectory 
simulations.  Errors  in  the  measurement  of  velocity  at  the  sampling  scale  were  considered. 
Dispersion  due  to  surface  water  currents  occurring  at  scales  smaller  than  the  sampling  scale 
(i.e.,  turbulence)  was  also  considered.  Methods  were  developed  to  quantify  the  error  and 
determine  dispersion  components  of  CODAR-derived  velocities.  Random  walk  (RW)  and 
random  flight  (RF)  Monte  Carlo  models  were  analyzed  and  compared.  SLDMB  trajectories 
were  then  modeled  incorporating  CODAR  and  STPS  velocity  data  using  both  the  RW  and  RF 
models. 

4.1.  Surface  current  dispersion 

The  motion  of  particles  (including  SAR  drift  objects)  is  caused  by  processes  that  vary  over  a 
wide  range  of  scales.  In  SAR  problems,  for  example,  processes  may  be  on  the  scale  of  the 
continental  shelf  or  down  to  the  scale  of  the  SAR  drift  object  itself.  Because  methods  are  not 
available  to  solve  for  motion  at  all  scales  simultaneously,  it  is  still  common  practice  (Taylor, 
1921)  to  describe  the  velocity  of  a  particle  in  terms  of  a  mean  field  component  and  a  turbulent 
component. 

Surface  current  mapping  radars  such  as  CODAR  can  provide  estimates  of  the  mean  field  velocity 
at  temporal  scales  of  1  hour  and  spatial  scales  of  approximately  1.5  km  and  6  km  for  standard- 
range  and  long-range  systems,  respectively.  This  radar-derived  mean  field  velocity  is  subject  to 
significant  uncertainties  and  errors,  which  effectively  constitute  a  separate  ‘error’  component  of 
the  total  velocity.  The  turbulent  component  of  velocity  incorporates  velocity  fluctuations  on  all 
scales  smaller  than  the  scale  of  the  CODAR-derived  mean  field  velocity.  For  a  detailed 
description  of  these  three  components  of  the  total  velocity,  see  (O'Donnell,  Ullman,  Edwards, 
Fake,  Allen,  2005)  and  Ullman  et  al,  (2005))  Appendices  C  and  D. 

For  reasons  of  practicality,  the  ‘error’  component  of  the  total  CODAR-derived  velocity  is 
combined  with  the  turbulent  component  of  the  velocity.  This  combination  is  partly  motivated  by 
the  fact  that  there  is  no  easy  way  to  separate  the  two  components  when  comparing  the  radar- 
based  mean  velocity  with  similar  estimates  from  drifters.  This  report  will  assume  that  this 
combined  error  plus  the  turbulent  velocity  component  can  still  be  described  by  models  of 
turbulence,  and  it  is  treated  as  such  in  this  study. 

The  dispersion  of  an  array  of  drifters  is  the  increase  in  mean  separation  of  neighboring  drifters 
over  time  due  to  the  effects  of  turbulence  and  spatial  gradients  (shear)  in  the  velocity  field.  To 
the  first  approximation,  dispersion  is  defined  in  terms  of  the  dispersion  coefficient,  K,  the 
velocity  variance,  a2,  and  the  autocorrelation  time  of  the  velocity  field,  x,  according  to  the 
following  equation  (Csanady,  1973): 

K=a\. 
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The  value  of  K  chosen  for  a  search  case  has  a  significant  impact  on  the  size  of  the  containment 
area.  As  an  example,  results  of  a  vertical  person  in  the  water  (PIW)  3-hour  drift  using  the  RW 
Monte  Carlo  method  and  dispersion  coefficients  of  5  m2/s  (the  value  originally  recommended  for 
the  operational  T&E)  and  140  m2/s  (Ullman’s  estimate  based  on  data  from  the  MAB  region) 
were  simulated  during  the  operational  test  and  evaluation  (T&E).  The  results  are  shown  in 
Figure  4-1.  The  containment  area  was  1.8  square  nmi  when  K  was  5  m2/s  and  34  square  nmi 
when  K  was  140  m  /s. 


Figure  4-1.  Containment  areas  resulting  from  different  dispersion  coefficients. 


In  practice,  dispersion  is  difficult  to  measure  directly.  Several  approaches  were  taken  to  quantify 
dispersion  in  the  BIS  and  the  MAB  CODAR  regions.  The  dispersion  coefficient  was  estimated 
directly  through  analysis  of  the  spread  of  individual  clusters  of  SLDMBs  with  time.  The 
dispersion  coefficient  can  also  be  calculated  as  the  product  of  the  velocity  variance  and  the 
velocity  autocorrelation  time.  Dispersion  coefficients  were  calculated  in  this  way  using  variance 
computed  for  CODAR  velocities  and  variance  computed  for  the  differences  between  CODAR 
and  SLDMB  velocities.  Dispersion  was  also  estimated  retrospectively  using  the  results  of  Monte 
Carlo  simulations  of  SLDMB  trajectories.  These  various  methods  resulted  in  estimates  of  the 
dispersion  coefficient  ranging  from  approximately  40  m2/s  to  700  m2/s  as  reported  in  (Spaulding, 
Isaji,  Hall,  Allen,  2005)  Appendix  E.  Ullman  et  al,  (2005)  used  the  variance  of  the  difference 
between  SLDMB  and  CODAR  velocities  (see  Appendix  B)  to  produce  estimates  of  K  for  the 
BIS  and  the  MAB  regions  that  are  presented  in  Table  4-1 .  This  approach  using  differences 
between  CODAR  and  SLDMB  velocities  is  considered  advantageous  conceptually  because  it 
most  closely  characterizes  the  sub-grid  scale  velocity  field  and  includes  the  measurement  error 
from  the  CODAR  system. 
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Table  4-1.  Summary  of  Dispersion  Coefficients  K  for  BIS  and  MAB  regions  calculated  using 
the  variance  of  the  difference  between  SLDMB  and  CODAR  velocities  in  the  BIS 


and  t 

re  MAB  regions. 

Region 

K-East-West  (m2/s) 

K-North-South  (m2/s) 

K-Average  (m2/s) 

BIS 

110 

80 

90 

MAB 

120 

160 

140 

4.2.  Monte  Carlo  characterization  of  dispersion 

The  random  walk  characterization  of  surface  current  dispersion  is  based  on  the  assumption  that 
the  non-tidal  residual  component  of  the  surface  current  vector  is  randomly  distributed.  This 
process  adds  random  noise  to  the  surface  current  vector  fields  to  reproduce  the  effects  of 
dispersion.  This  method  is  easily  implemented  and  has  seen  frequent  use  in  past  Monte  Carlo 
search  planning  tools.  In  the  Random  Flight  characterization  of  surface  currents,  dispersion  is 
modeled  by  drawing  values  that  are  random,  yet  are  partially  dependent  on  the  previous  value  of 
the  surface  current  vector.  Each  new  vector  calculated  by  the  RF  approach  is  therefore  partially 
correlated  with  the  previous  value.  As  new  vectors  are  calculated  on  succeeding  time  steps, 
positions  at  later  time  steps  become  less  dependent  on  earlier  values.  The  system  “decorrelates”, 
approximately  at  an  exponential  rate  with  the  number  of  time  steps.  The  characteristic  time 
corresponding  to  that  exponential  decay  rate,  is  referred  to  as  the  “e-folding”  time  or 
autocorrelation  time. 

In  practice,  the  RW  model  treats  dispersion  through  the  dispersion  coefficient,  K.  The  RW 
approach  handles  dispersion  (S2),  or  the  mean  squared  separation  of  particles,  as  linear  growth 
with  time  (t)  at  the  rate  S2  =  2Kt.  With  the  RF  model,  for  a  given  value  of  dispersion,  the  mean 
squared  particle  separation  grows  at  a  slower  rate  of  time,  S2  =  2Kt  -  2KT(l-e't/T),  where  T  is  the 
characteristic  autocorrelation  time  of  the  velocity  field  (Spaulding  et  al,  2005).  For  a  given  K, 
the  dispersion  of  a  particle  field  predicted  by  the  RW  approach  will  be  larger  than  one  predicted 
by  the  RF  approach.  For  detailed  methodologies  of  these  models,  also  see  Griffa  (1996)  or 
Berloff  and  McWilliams  (2002). 

A  comparison  of  analytical  solutions  for  the  RW  and  RF  models  shows  that  for  a  given  value 
of  K,  the  RW  model  results  in  more  dispersion  (and  a  larger  predicted  search  area)  than  the  RF 
model  (Figure  4-2).  This  difference  is  particularly  pronounced  at  times  shorter  than  the  velocity 
autocorrelation  time  (T=3  hr).  During  this  period,  the  size  of  the  search  area  predicted  by  the 
RW  model  increases  more  rapidly  than  the  size  predicted  by  the  RF  model.  At  times  greater  than 
the  velocity  autocorrelation  time,  the  search  areas  predicted  by  the  two  models  increase  at  the 
same  rate. 
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Time  (h  ouix) 

Figure  4-2.  Analytical  solution  of  predicted  search  areas  (separation)  as  a  function  of 

time  for  the  RW  and  RF  models.  The  difference  in  search  areas  predicted  by 
the  two  models  is  small  relative  to  the  overall  range  of  potential  search  areas 
due  to  uncertainty  in  K.  For  example,  RW  and  RF  separation  initially 
increases,  but  levels  off  at  a  constant  value  for  a  given  K. 

While  the  RF  approach  is  conceptually  better,  the  difference  between  the  search  areas  obtained 
from  RW  and  RF  models  is  small  relative  to  the  difference  resulting  from  uncertainty  in 
dispersion  ( K ).  The  value  selected  for  K  emerges  as  the  critical  factor  in  defining  the  size  and 
growth  of  the  search  area  in  a  SAR  case.  Figure  4-2  therefore  underscores  the  significance  of 
having  a  representative  value  of  K  for  SAR  planning.  If  K  is  overestimated,  the  search  area  will 
be  too  large,  resulting  in  a  waste  of  search  effort.  If  K  is  underestimated,  the  search  area  will  be 
too  small,  increasing  the  probability  that  the  search  object  will  be  outside  of  the  search  area  and 
increasing  the  potential  that  the  search  object  will  not  be  found. 

4.3.  Monte  Carlo  simulation  results 

Simulations  were  conducted  to  evaluate  the  ability  of  the  two  Monte  Carlo  models  to  predict 
search  areas  using  CODAR  and  STPS  velocity  data  in  conjunction  with  the  dispersion 
parameters  obtained  from  this  study  (using  the  East-West  and  North-South  dispersion 
coefficients  for  BIS  and  MAB  in  Table  4-1).  Search  areas  were  predicted  for  the  individual 
24-hour  BIS  SLDMB  trajectories.  A  total  of  159  independent  drifter  trajectories  were  modeled 
in  the  MAB;  72  drifter  trajectories  were  modeled  in  the  BIS.  Figure  4-3  shows  an  example  of  an 
actual  SLDMB  trajectory  with  a  predicted  search  area  derived  from  pseudo-drifter  trajectories 
derived  from  a  RF  simulation.  The  95  percent  containment  region  (predicted  search  area)  was 
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estimated  by  computing  a  two-dimensional  histogram  of  the  pseudo-drifter  locations,  sorting  the 
spatial  bins  based  on  the  number  of  pseudo-drifters  contained,  and  then  summing  over  the  bins 
until  95  percent  of  the  total  was  obtained.  This  approach  allowed  the  containment  area  to  be 
more  closely  tailored  to  the  distribution  of  the  modeled  pseudo-drifters  than  a  simple  circular 
area  would  be.  The  performance  of  the  RW  and  RF  approaches  was  compared  by  calculating  the 
percent  of  cases  when  the  actual  SLDMB  was  located  in  the  predicted  search  area  between  1  and 
24  hours.  If  the  SLDMB  was  located  within  the  predicted  search  area  (defined  as  the  95  percent 
containment  region),  the  prediction  was  considered  a  success. 


Drifter  36446. 12/1 7/2Q02 1 1 :00 


Longilude  (’W) 

Figure  4-3.  Illustration  showing  actual  24  hr  SLDMB  track  (red)  and  predicted  drifter 
trajectory  (green).  Small  blue  dots  represent  the  positions  of  the  individual 
particles  used  in  a  series  of  Random  Flight  Monte  Carlo  simulations.  The 
95  percent  containment  region  is  shown  as  gray-shaded  rectangles. 


Using  the  drifter  trajectories  described  above  and  corresponding  CODAR  data,  Ullman  et  al, 
(2005)  calculated  dispersion  coefficients  for  RF  and  RW  trajectory  simulations  in  BIS  and  MAB. 
The  dispersion  coefficients  were  based  on  the  variance  and  autocorrelation  times  of  the  nontidal 
velocity  plus  measurement  error  time  series.  Ullman  et  al,  (2005)  found  that  the  area  containing 
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95  percent  of  SLDMB  trajectories  predicted  by  the  RF  approach  contained  drifters  with  a 
probability  of  about  90  percent.  On  the  other  hand,  a  much  smaller  95  percent  containment  area 
was  predicted  by  the  RW  approach  that  contained  the  drifter  in  only  60  percent  to  80  percent  of 
the  cases.  The  disparity  was  attributed  to  an  underestimate  of  the  autocorrelation  time  estimate 
used  for  K  in  the  RW  calculations.  When  K  was  changed  to  a  value  consistent  with  the  one  used 
for  the  RF  simulations,  the  two  approaches  provided  consistent  results  (Figure  4-4),  particularly 
at  times  greater  than  the  velocity  autocorrelation  time.  This  behavior  was  demonstrated  in  the 
analytic  solutions  for  the  two  models  (Figure  4-2). 


Time  (h) 

Figure  4-4.  A  comparison  of  RW  and  RF  Monte  Carlo  simulation  results  for  the  MAB  region. 

For  each  hourly  prediction  time,  the  percent  of  cases  where  the  actual  drifter  location 
at  that  time  fell  within  the  estimated  95  percent  confidence  region  is  plotted. 

4.4.  Discussion  of  Monte  Carlo  modeling  results 

The  Monte  Carlo  exercises  were  somewhat  inconclusive  in  that  there  was  some  disagreement 
over  the  relative  merits  of  the  Random  Walk  and  Random  Flight  approaches.  For  short  Search 
Action  Plan  drifts  (i.e.,  on  the  order  of  the  autocorrelation  time),  the  RW  model  produced  a 
larger  predicted  search  area  than  the  RF  model  and  consequently  was  more  likely  to  contain  the 
target.  As  drift  times  lengthen,  the  relative  difference  in  predicted  search  area  between  the  two 
models  decreased,  and  the  performance  of  the  two  models  was  approximately  equivalent. 
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The  overarching  issue,  however,  appears  to  be  that  the  value  selected  for  K  exerts  a  major 
influence  on  the  outcome  of  the  initial  search  area  that  would  be  estimated.  Values  of  K  between 
40  to  700  m  /s  produced  separation  estimates  from  7  to  120  square  kilometers  in  24  hours,  based 
upon  Spaulding  et  al,  (2005)  (Appendix  E),  equation  4  (RW  Separation  area  =  2  x  K  x  time). 
Using  the  difference  between  SLDMB  and  CODAR  velocities  in  the  BIS  region,  K  was 
estimated  to  be  90  m2/s.  In  the  MAB  region,  K  was  estimated  to  be  140  m2/s.  Thus,  the 
separation  area  after  24  hours  of  drift  for  the  BIS  area  is  expected  to  be  16  square  kilometers  and 
24  square  kilometers  for  the  MAB  region. 

The  most  practical  ways  to  estimate  K  are  determining  the  dispersion  parameters  (i.e.,  velocity 
variance  and  autocorrelation  time)  from  CODAR  velocity  data,  SLDMB  velocity  data,  or  from 
the  difference  between  CODAR  and  SLDMB  velocities.  Using  the  difference  between  CODAR 
and  SLDMB  velocities  most  closely  characterizes  the  sub-grid  scale  velocity  field  and  includes 
the  measurement  error  from  the  CODAR  system.  While  analyzing  the  dispersion  of  the 
individual  SLDMB  clusters  is  the  most  direct  way  to  measure  dispersion,  this  process  is 
problematic  in  that  large  numbers  of  SLDMBs  (i.e.,  on  the  order  of  10)  need  to  be  deployed  and 
allowed  to  drift  for  a  significant  period  of  time  in  order  for  the  resulting  statistics  to  be  truly 
robust.  Completing  this  task  for  the  CG’s  area  of  responsibility  would  be  a  large  task. 

Dispersion  parameters  vary  both  spatially  and  temporally,  and  the  degree  to  which  dispersion 
parameters  can  be  usefully  differentiated  depends  on  both  the  spatial  and  temporal  density  of 
data  used  to  estimate  them.  Improved  estimates  of  dispersion,  including  the  possibility  of  some 
spatial  refinement,  can  be  achieved  incrementally  through  additional  SLDMB  deployments.  A 
system  could  be  established  to  catalog  SLDMB  and  CODAR  data  for  individual  CODAR 
regions.  By  maintaining  a  database  with  all  SLDMB  trajectories  through  a  given  CODAR 
region,  along  with  concurrent  CODAR  velocities,  estimates  could  be  refined  as  necessary  in  the 
future.  Calculations  could  be  performed  on  a  periodic  basis  (e.g.,  annually).  As  new  data  are 
added  to  the  calculations,  the  accuracy  of  the  estimates  of  K  would  approach  their  representative 
mean  values  over  a  period  of  several  years. 
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5.  EXTENSION  OF  THE  PROTOTYPE  SAROPS  ENVIRONMENTAL 
DATA  RETRIEVAL  SUBSYSTEM  TO  EXTRACT  STANDARD  CODAR 
DATA  FIELDS  AND  THEIR  ASSOCIATED  STPS  FORECASTS 

5.1.  Prototype  SAROPS  EDS 

SAROPS  is  the  new  CG  SAR  system  currently  being  developed  by  Northrop  Grumman,  Applied 
Science  Associates  (ASA),  and  Metron.  SAROPS  contains  three  main  components: 
ArcGIS/Commercial  Joint  Mapping  Tool  Kit  (CJMTK)-based  graphical  user  interface  (GUI), 
Environmental  Data  Server  (EDS),  and  Simulator  engine  that  performs  the  particle  motion  and 
search  optimization. 

At  the  time  of  this  project,  SAROPS  was  still  under  development;  the  alpha  version  of  SAROPS 
was  delivered  in  March  2005,  and  final  release  is  not  expected  until  2006.  As  a  replacement  for 
the  SAROPS  application,  a  modified  version  of  SARMAP,  a  commercially  available  SAR 
application  used  by  international  CG  agencies,  was  substituted.  SARMAP  includes  geographic 
information  system  (GlS)-based  tools  for  developing  SAR  cases,  gathering  and  entering 
environmental  data,  running  International  Aeronautical  and  Maritime  Search  and  Rescue 
(IAMSAR)  or  Monte  Carlo  simulations,  visualizing  SAR  results,  and  adding  Search  and  Rescue 
Unit  (SRU)  patterns. 

Even  though  SAROPS  was  not  available,  a  preliminary  version  of  the  SAROPS  EDS  was 
available  for  evaluation.  SARMAP  was  customized  to  connect  to  the  EDS  to  access 
environmental  data.  This  connection  enabled  evaluation  of  the  SAROPS  client  communications 
with  the  EDS  by  using  the  modified  SARMAP  software  as  a  proxy  for  SAROPS. 

Primary  data  served  by  the  EDS  are  current  and  wind  data.  Additional  environmental  parameters 
are  expected  to  be  available  in  the  future.  Environmental  data  are  based  on  the  SAR  case 
information.  SAROPS  determines  the  area  of  interest  (AOI)  and  time  of  interest  before  making  a 
request  to  the  EDS  for  data  (currents  and  winds)  that  meet  the  AOI  and  time  of  interest. 

The  initial  developmental  prototype  EDS  version  retrieved,  archived,  and  distributed  data  from 
NOAA  PORTS  sites,  NOAA  CO-OPS,  selected  NDBC  stations,  and  selected  New  England 
USGS  river  gauges.  As  part  of  the  continuing  SAROPS  development  effort,  the  EDS  was 
modified  to  add  CODAR  and  SLDMB  data.  The  EDS  is  made  up  of  two  main  components:  a 
catalog  server  and  an  EDS  Web  server. 

5.1.1.  Catalog  server 

The  catalog  component  of  the  EDS  contains  the  data  sources  available  to  the  EDS.  Both 
HyperText  transfer  protocol  (HTTP)  and  file  transfer  protocol  (FTP)  were  required  to  retrieve 
information  from  the  data  sources.  The  significant  information  contained  in  the  catalog  are: 

•  Geographic  extent  of  the  data  source, 

•  Time  window  of  data  available  for  a  data  source,  and 

•  Structure  and  format  of  the  source  files. 
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On  a  routine  basis,  based  on  source  data  update  frequency,  the  catalog  server  retrieves  data  from 
different  sources  and  makes  these  data  available  to  the  SAROPS  user.  As  new  data  are  collected 
and  stored,  the  catalog  is  updated  to  reflect  the  most  recent  “time-stamp”  of  available  data. 
Archives  are  kept  so  that  SAROPS  users  can  retrieve  up  to  six  months  of  historical  data  for  case 
reconstruction. 

Environmental  data  are  highly  heterogeneous;  data  providers  supply  their  data  in  a  variety  of 
formats.  In  some  cases,  the  catalog  server  converts  the  incoming  data  to  make  access  and 
archiving  simpler  and  to  provide  compatibility  with  common  units,  coordinates,  etc. 

The  modified  EDS  collects  and  distribute  four  data  products: 

•  UCONN  CODAR  and  STPS  data 

The  EDS  catalog  retrieved  the  BIS  CODAR  and  STPS  files  on  an  hourly  basis  from  the 
UCONN  server  using  FTP  and  stored  files  on  the  EDS. 

•  Rutgers  CODAR  and  STPS  data 

The  EDS  catalog  retrieved  the  MAB  CODAR  files  on  an  hourly  basis  from  the  Rutgers 
server  using  FTP  and  stored  files  on  the  EDS.  The  MAB  STPS  files  were  retrieved  on  an 
hourly  basis  from  the  UCONN  server-using  FTP,  and  files  were  stored  on  the  EDS. 

•  Selected  northeast  NDBC  wind  station  data 

The  catalog  stored  uniform  resource  locator  (URL)  information  for  selected  NDBC 
stations  (Table  5-1)  and  collected  the  data  using  HTTP  on  an  hourly  basis,  parsed  the 
HyperText  markup  language  (HTML)  data,  and  appended  data  for  each  station  into  a 
delimited  text  file. 

Table  5-1.  NDBC  stations  accessed  by  prototype  EDS. 

NDBC  Station  List 

44004  -  HOTEL 
44005 -GULF  OF  MAINE 

44007  -  PORTLAND,  ME 

44008  -  NANTUCKET 

44009  -  DELAWARE  BAY _ 

4401 1  -  George’s  Bank 

44013 -BOSTON,  MA _ 

44017  -  Montauk  Point,  NY 

4401 8 -SE  Cape  Cod _ 

44025  -  LONG  ISLAND _ 

44027  -  Jonesport,  ME _ 

ABAN6  -  Alexandria  Bay,  NY _ 

ALSN6  -  Ambrose  Light,  NY _ 

BUZM3  -  Buzzard’s  Bay,  MA 

IOSN3  -  Isle  of  Shoals,  NH _ 

MDRM1  -  Mt  Desert  Rock,  ME _ 

MISM1  -  Matinicus  Rock,  ME 
SUPN6  -  Superior  Shoals,  NY 

THIN6  -  Thousand  I  Bridge,  NY _ 

TPLM2  -  Thomas  Point,  MD _ 

•  SLDMB  tracks  (from  Operations  System  Center  (OSC)  Martinsburg). 

The  catalog  retrieved  the  latest  SLDMB  tracks  from  an  American  Standard  Code  for 
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Information  Exchange  (ASCII)  comma-delimited  text  file  stored  on  an  FTP  site  at 
OSC  Martinsburg. 

The  EDS  version  implemented  for  this  project  is  the  first  draft  (alpha)  of  the  EDS  that  will  be 
installed  for  the  SAROPS  project.  Evaluation  of  the  EDS  shows  that  the  use  of  Network 
Common  Data  Form  (NetCDF)  as  a  data  transport  format  is  quite  efficient  and  also  allows  the 
EDS  to  transition  to  use  of  Open-source  Project  for  a  Network  Data  Access  Protocol 
(OPeNDAP)  in  the  future.  There  are  many  issues  involved  in  the  long-term  evolution  and 
improvement  of  the  EDS  including: 

•  Server  side  filtering  (spatial  and  temporal)  based  on  user  data  requests; 

•  Server  side  aggregation  to  combine  data  across  space  and  time; 

•  Server  side  quality  assurance/quality  control  (QA/QC)  processes;  and 

•  Server  side  data  blending,  objective  analysis,  and  data  assimilation  techniques. 

5.1.2.  EDS  Web  Services 

The  EDS  Web  services  (.NET  based)  wait  for  requests  from  the  SAROPS  client.  The  requests 
come  in  extensible  Markup  Language  (XML)  format  through  Simple  Object  Access  Protocol 
(SOAP).  When  a  request  is  received  for  meteorological  and  oceanographic  data  for  a  particular 
area  of  interest  (AOI)  and  time  of  interest,  the  Web  service  processes  the  appropriate  data  from 
the  catalog  and  returns  data  to  the  SAR  client.  A  schematic  of  this  process  is  included  in  Figure 
5-1.  The  modified  EDS  Web  services  manage  the  following:  CODAR/STPS,  NDBC  winds,  and 
SLDMB  data. 


Figure  5-1.  Schematic  of  SARMAP  (SAROPS)  clients  connecting  to  EDS  for  environmental 
data  needs. 
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5.2.  CODAR/STPS 


A  user  requesting  water  current  data  from  the  EDS  can  select  from  one  of  three  choices  as  shown 
in  Figure  5-2. 

•  CODAR 

This  process  gives  the  user  the  “best”  available  data  by  combining  the  latest  available 
CODAR  (observation)  data  with  STPS  data  for  a  forecast  period.  Although  these  data 
are  uniform,  the  operational  test  and  evaluation  (T&E)  users  noted  that  this  solution  did 
not  allow  them  to  know  whether  a  particular  time  step  was  represented  by  CODAR 
(actual  observation)  or  STPS  (forecast).  As  a  result,  the  users  were  uncertain  as  to  how 
much  confidence  to  put  in  the  results. 

•  CODAR  RAW 

This  process  gives  the  user  only  CODAR  data  and  does  not  fill  any  gaps  with  STPS  or 
provide  STPS  forecast  data. 

•  STPS 

This  process  gives  the  user  only  STPS  data;  no  actual  CODAR  observation  data  are  used. 


Figure  5-2.  Example  SARMAP  Graphical  User  Interface  (GUI)  connecting  to  EDS  for  user 
selection  of  environmental  data. 

5.3.  NDBC  winds 

In  response  to  a  user  request  for  wind  data,  the  EDS  provides  a  list  of  NDBC  wind  stations  in  the 
region  (Figure  5-3).  The  user  selects  a  wind  station,  and  the  wind  time  series  for  a  SAR  case  is 
automatically  loaded  from  the  NDBC  data. 


5-4 


Figure  5-3.  Example  SARMAP  Graphical  User  Interface  (GUI)  connecting  to  the  EDS  for  user 
selection  of  NDBC  wind  data. 

The  NDBC  data  are  stored  on  the  EDS  as  Coastmap  Text  Format  (CTF)  files,  are  parsed  for  the 
time  of  interest,  and  then  returned  to  the  client  as  a  binary  wind  fde  (.WNE).  The  CTF  file 
format  is  a  flexible  space-delimited  test  format  useful  for  storing  time  series  data  (e.g.,  wind 
speed  and  direction,  tide  height,  water  temperature)  for  a  single  geographic  location.  Metadata 
are  included  in  the  file  header,  and  all  data  at  a  given  time  are  listed  on  single  lines  with  a  date 
and  time  stamp. 

5.4.  SLDMB  data 

In  response  to  a  user  request  for  SLDMB  data,  the  EDS  processes  the  ASCII  file  containing  all 
of  the  SLDMB  locations  in  near  real  time  and  returns  an  Environmental  Systems  Research 
Institute,  Inc.  (ESRI)  shape  (SHP)  file  containing  the  buoy  tracks  for  the  time  of  interest.  This 
SHP  file  is  then  displayed  on  the  SARMAP  chart/map. 

5.5.  Conclusions  on  extension  of  the  prototype  SAROPS  environmental  data 
retrieval  subsystem 

This  task  demonstrated  that  it  will  be  possible  for  the  SAROPS  EDS  to  collect  and  distribute 
CODAR  current  data  and  STPS  current  predictions  from  UCONN  and  Rutgers,  wind 
observations  from  NDBC  buoy  sites,  and  SLDMB  data  from  OSC  Martinsburg.  It  will  also  be 
possible  for  remote  users  to  request  and  download  these  data  in  near  real-time  for  use  in  Search 
Action  Plans  (SAPs). 
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6.  OPERATIONAL  TEST  AND  EVALUATION 


CODAR  and  STPS  forecasts  offer  greater  temporal  and  spatial  resolution  than  those  used  in 
C2PC/C ASP/joint  automated  worksheets  (JAWS)  software  and  may  be  useful  in  developing 
SAPs  that  confine  searches  to  smaller  areas  and  increase  the  probability  of  success  (POS). 

During  the  operational  T&E,  controllers  had  the  opportunity  to  compare  SAP  data  developed 
from  CODAR/STPS  current  data  with  comparable  data  from  C2PC/JAWS.  They  also  had 
opportunities  to  provide  suggestions  for  improvements  to  CODAR/STPS  products  and  the 
software  that  accessed  these  products. 

The  two  objectives  of  the  30-day  operational  T&E  were: 

1)  Validate  the  utility  of  using  near-real-time  CODAR  and  STPS  forecast  data  interpolated  to  a 
small  geographic  area  to  develop  a  search  action  plan,  and 

2)  Assess  the  level  of  confidence  that  CG  operations  personnel  (i.e.,  controllers)  have  in  using 
CODAR/STPS  data  to  develop  a  search  action  plan. 

6.1.  Test  description 

CG  Group  Moriches  (GRUMOR)  hosted  the  operational  T&E.  Operations  Center  (OPCEN) 
controllers  were  trained  on  the  use  of  the  modified  SARMAP  software  and  its  wind  and  data 
retrieval  functions  to  access  the  EDS  where  CODAR/STPS  data  are  posted.  All  Group 
controllers  participated  during  the  30-day  evaluation  period  by  accessing  CODAR/STPS  data  to 
execute  simulated  test  cases,  complete  data  forms,  and  share  their  comments  and  suggestions 
during  interviews  with  the  test  observer. 

Pretest  training  was  provided  from  20-3 1  October  2004  on  an  ad  hoc  basis  at  times  convenient  to 
the  controllers.  Each  controller  received  approximately  four  hours  of  one-on-one  instruction 
during  the  initial  week  of  training.  During  the  first  week  of  the  operational  T&E,  additional 
follow-on  training  was  conducted  to  ensure  that  the  controllers  were  using  the  software  correctly, 
had  no  questions,  and  understood  their  roles  during  the  test  period. 

The  operational  T&E  was  conducted  from  1-29  November  2004.  During  this  period,  four 
separate  predefined  six-hour  observation/interview  sessions  were  conducted  by  a  test 
observer/data  collector. 

During  the  operational  T&E,  the  following  three  SLDMBs  were  deployed  as  illustrated  in 
Figure  6-1  to  provide  simulated  zero-leeway  search  objects. 

1)  #38968  on  27  October  2004  at  041°  07.00'  N,  071°  42.00'  W. 

2)  #38708  on  3  November  2004  at  040°  37.90'  N,  072°  28.90'  W. 

3)  #39001  on  18  November  2004  at  040°  26.23'  N,  073°  18.18'  W. 
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SLDMB  #38968  stayed  in  the  BIS  CODAR  region  from  27-31  October,  drifted  south  then  east 
during  a  period  of  high  wind,  and  continued  southwest  until  it  entered  the  MAB  region  on 
15  November.  SLDMB  #38968  then  remained  in  the  MAB  region  for  the  duration  of  the  test. 
SLDMB  #38708  remained  on  the  northeast  edge  of  the  MAB  region  from  3-5  November,  then 
drifted  east  during  a  period  of  high  wind,  continued  southeast  until  it  re-entered  the  MAB  region 
on  14  November,  and  then  remained  in  the  MAB  region  for  the  duration  of  the  test. 

SLDMB  #39001  remained  on  the  northern  border  of  the  MAB  region  from  its  deployment  on 
18  November  throughout  the  operational  T&E  period. 

GRUMOR  OPCEN  was  contacted  once  per  shift  and  supplied  with  simulated  case  data  that 
could  be  used  to  exercise  the  modified  SARMAP  software  and  the  C2PC/JAWS  software.  The 
data  were  developed  from  the  SLDMB  data  and  specified  the  target,  the  Greenwich  Mean  Time 
(GMT)  of  the  LKP  to  be  used,  the  latitude/longitude  of  the  LKP,  and  the  length  of  time  of  the 
drift  for  the  test  case. 

All  watchstanders  were  asked  to  run  the  simulated  case  twice,  once  using  modified  SARMAP 
and  once  using  C2PC/JAWS.  They  were  requested  to  archive  all  computer-generated  files 
required  by  the  modified  SARMAP  software  so  that  the  cases  could  be  reconstructed  following 
the  test,  and  to  save  case  summaries  in  text  file  format  from  both  the  C2PC/JAWS  and  the 
modified  SARMAP  software.  After  running  a  simulated  SAR  case,  and  waiting  for  the  SLDMB 
trajectory  for  the  time  period  of  interest  to  be  completed,  the  controller  could  access  the  SLDMB 
data  to  determine  the  accuracy  of  the  modeled  containment  area. 
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In  addition,  controllers  were  also  asked  to  use  SARMAP  to  develop  a  search  action  plan  for  any 
actual  offshore  SAR  cases  that  developed  during  their  watch.  These  cases  could  be  run 
whenever  it  was  convenient.  As  it  turned  out,  there  were  no  cases  appropriate  to  this  purpose 
during  the  operational  T&E  period. 

6.2.  Controller  feedback 

In  general,  the  GRUMOR  controllers  were  positive  and  enthusiastic  about  accessing 
CODAR/STPS  data  using  the  modified  SARMAP  software.  Collected  feedback  is  presented  in 
Table  6-1. 


Table  6-1.  Controller  feedback. 


Feedback 

Applicability 

Feedback 

Modified  SARMAP 

Simple  to  use. 

Requires  fewer  clicks  of  the  mouse;  clickable  functions  save  time. 

There  are  fewer  data  to  input. 

CODAR/STPS  data 

Not  appropriate  for  surf-zone  SAR  cases. 

MAB  and  BIS  regions  do  not  provide  sufficient  coverage  for  the  GRUMOR  area 
of  responsibility  (AOR). 

CODAR  data  encompassed  the  target  intermittently  because  the  size  and  shape 
of  the  CODAR  coverage  area  varied  with  time. 

It  is  impossible  to  tell  if  the  displayed  current  data  are  from  CODAR  or  from 

STPS. 

Wind  data 

When  predicted  wind  is  applied  to  STPS  forecast  current  data,  it  is  also 
inappropriately  applied  to  CODAR  currents. 

Current  data 

Cannot  mix  CODAR  analysis  with  tidal  predictions  or  datum  marker  buoy  (DMB) 
data  with  STPS  forecasts. 

Search  area  results 

The  SARMAP  containment  area  is  small  compared  with  the  C2PC/JAWS  Total 
Probable  Error  (TPE). 

The  target  was  not  always  contained  in  the  containment  area. 

Display 

The  ability  to  view  the  object  drift  on  the  map  and  display  the  surrounding 
current  vectors  are  good  features. 

Modified 

SARMAP/EDS 

interaction 

When  the  incident  site  is  outside  of  a  CODAR  analysis  coverage  area,  the 
software  should  let  the  user  know  that  the  data  may  not  be  appropriate. 

Modified  SARMAP  uses  the  four  nearest  vectors  (regardless  of  how  far  away 
they  are)  to  define  the  current  field. 

6.3.  SAROPS  and  EDS  interaction 

Decision  automation  was  an  item  of  particular  interest.  The  controllers  recommend  that  the  SAP 
software  should  suggest  an  NDBC  station,  a  CODAR  analysis  area,  and  parameter  values  (such 
as  dispersion,  safety,  and  drift)  based  on  the  incident  location  to  be  used  as  input  to  the  SAP 
processing.  At  the  same  time,  they  pointed  out  that  C2PC/JAWS  currently  allows  them  to  select 
and  view  the  NDBC  buoy  to  ensure  that  it  is  operational  and  that  the  data  from  the  buoy  are  in 
the  expected  range  before  the  data  are  imported  to  the  SAR  software. 
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6.4.  Test  case  comparisons 

6.4.1.  STPS  forecast  versus  CODAR  analysis 

The  procedures  described  in  Section  6.1  yielded  47  test  cases  for  evaluation.  All  test  cases  were 
later  reviewed  to  compare  separation  distances  to  the  SLDMB  drifter  when  the  different  models 
were  used  and  to  associate  the  quality  of  the  modified  SARMAP  results  with  the  proximity  to 
valid  CODAR/STPS  data.  Watchstander  workloads  during  the  tests  precluded  an  accurate 
comparison  of  the  two  models.  Time  expended  accessing  the  SLDMB  data,  developing  case 
parameters,  reporting  the  case  to  the  OPCEN,  and  executing  the  case  at  an  opportune  time  often 
caused  cases  to  be  run  later  in  or  after  the  STPS  prediction  period.  Once  the  prediction  time 
became  a  CODAR  observation,  STPS  values  for  that  time  were  overwritten  as  CODAR  data  and 
a  new  set  of  STPS  values  were  created.  When  this  occurred,  these  cases  were  based  only  on 
CODAR  data  and  not  on  any  STPS  predictions,  preventing  an  accurate  comparison  of  the  two 
approaches. 

Some  cases  were  based  on  time  periods  when  the  target  SLDMB  was  not  in  the  CODAR/STPS 
data  coverage  area.  For  most  of  the  cases,  the  target  SLDMB  was  in  the  area,  but  the  CODAR 
data  coverage  area  only  intermittently  included  it.  This  was  because  returns  measured  from  the 
HF  radar  sensors  propagate  in  the  atmosphere  inconsistently,  especially  at  night  and  cause  the 
reported  current  field  size  to  vary  by  the  hour.  The  STPS  forecast  data,  based  on  averaged 
CODAR  data,  has  a  larger  area  of  consistent  coverage. 

6.4.2.  CODAR/STPS  versus  C2PC/JAWS  comparisons  inside  the  CODAR/STPS  coverage 
area 

Eleven  simulated  SAR  cases  were  completely  inside  the  CODAR  analysis/STPS  forecast  area. 
The  International  Aeronautical  and  Maritime  Search  and  Rescue  (IAMSAR)  model  was  used  to 
match  the  C2PC/JAWS  model  currently  used  at  GRUMOR.  Note:  One  case  is  included  twice 
(once  for  STPS,  and  once  for  CODAR).  Modified  SARMAP  predictions  from  these  cases  were 
compared  with  C2PC/JAWS  predictions.  Results  are  shown  in  Figure  6-2.  Statistical 
parameters  are  provided  in  Table  6-2.  The  average  ranges  (separation  distances)  from  the 
centers  of  the  left  and  right  IAMSAR-modeled  position  estimates  to  the  SLDMB  positions  at  the 
ends  of  the  drifts  was  1.3  nmi  when  the  modified  SARMAP  with  CODAR/STPS  data  were  used 
compared  with  2.8  nmi  when  C2PC/JAWS  was  used.  More  variability  was  also  evident  using 
the  C2PC/JAWS  model:  1 1.0  nmi  compared  with  1.4  nmi  (SARMAP  value). 
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Drift  Duration  (hrs) 

Figure  6-2.  Comparison  of  CODAR/STPS  and  C2PC/JAWS  results  inside  the  CODAR/STPS 
coverage  area. 


Table  6-2.  Statistical  comparison  of  CODAR/STPS  versus  C2PC/JAWS  inside  the 
_ CODAR/STPS  coverage  area. _ 


Model 

Mean  Range  to  Drifter  (nmi) 

Standard  Deviation  (nmi) 

CODAR/STPS 

1.3 

1.4 

C2PC/JAWS 

2.8 

11.0 

6.4.3.  CODAR/STPS  versus  C2PC/JAWS  comparisons  outside  the  CODAR/STPS 
coverage  area 

Fourteen  simulated  SAR  cases  were  completely  outside  the  CODAR  analysis/STPS  forecast 
area.  The  IAMSAR  model  was  used  to  match  the  C2PC/JAWS  model  currently  used  at 
GRUMOR.  When  the  drifter  was  outside  the  CODAR/STPS  coverage  area,  the  modified 
SARMAP  software  extrapolated  from  the  nearest  available  CODAR/STPS  data.  Trajectories 
predicted  outside  of  the  CODAR/STPS  coverage  are  expected  to  be  less  accurate  because  the 
current  field  is  not  spatially  uniform  and  because  the  CODAR  data  used  in  the  extrapolation 
most  likely  had  higher  associated  uncertainties  as  they  are  on  the  fringe  of  the  coverage  area. 
Modified  SARMAP  predictions  from  these  cases  were  compared  with  C2PC/JAWS  predictions. 
Results  are  shown  in  Figure  6-3.  Statistical  parameters  are  provided  in  Table  6-3.  The  average 
range  from  the  centers  of  the  left  and  right  IAMSAR-modeled  position  estimates  to  the  SLDMB 
positions  at  the  ends  of  the  drifts  was  2.4  nmi  when  the  modified  SARMAP  with  CODAR/STPS 
data  were  used  compared  with  2.6  nmi  when  C2PC/JAWS  was  used.  More  variability  was 
evident  using  the  CODAR/STPS  data,  2.8  nmi  compared  with  0.7  nmi. 


6-5 


Drift  Duration  (hrs) 

Figure  6-3.  Comparison  of  CODAR/STPS  and  C2PC/JAWS  results  outside  the  CODAR/STPS 
coverage  area. 


Table  6-3.  Statistical  comparison  of  CODAR/STPS  versus  C2PC/JAWS  outside  the 


CODAR/STPS  coverage  area. 


Model 

Mean  Range  to  Drifter  (nmi) 

Standard  Deviation  (nmi) 

CODAR/STPS 

2.4 

2.8 

C2PC/JAWS 

2.6 

0.7 

6.5.  Summary  of  the  CODAR/STPS  operational  T&E 

•  The  BIS  and  the  MAB  CODAR  sites  were  fully  functional  during  the  entire  operational 
T&E  period.  This  demonstrated  the  operational  reliability  of  these  automated  processing 
systems. 

•  CG  GRUMOR  controllers  demonstrated  interest  and  confidence  in  using  CODAR 
surface  currents  measured  in  near  real  time  and  STPS  estimated  currents  developed  from 
CODAR  data  to  develop  SAPs. 

•  Use  of  the  modified  SARMAP  to  interface  with  the  EDS  raised  issues  that  should  be 
considered  in  the  development  of  SAROPS,  e.g.,  decision  automation. 

•  On  the  basis  of  a  very  limited  data  set,  the  modified  SARMAP  model  using  CODAR  and 
STPS  data  predicted  target  drift  of  simulated  PIW  targets  more  accurately  than  the 
C2PC/JAWS  model.  Caution  should  be  used  when  the  drifter  is  outside  the 
CODAR/STPS  coverage  area  as  the  results  may  be  based  on  incorrect  current  data  both 
because  the  current  may  be  very  different  at  the  drift  location  than  it  is  within  the 
CODAR/STPS  coverage  area  and  also  because  the  extrapolated  data  probably  have 
higher  associated  uncertainties  than  more  included  data. 
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7.  SUMMARY  AND  RECOMMENDATIONS 


7.1.  STPS  modifications 

The  STPS,  which  had  been  developed  for  the  BIS  where  surface  currents  were  considered  to  be 
tidally  dominated,  was  adapted  to  the  MAB,  where  local  winds  were  considered  to  dominate 
circulation.  Several  enhancements  were  made  to  the  STPS.  The  calculations  were  modified  to 
factor  in  the  localized  wind-current  relationships  at  each  grid  point,  instead  of  a  region- wide 
average.  The  cross-correlations  between  the  east  and  north  velocity  components  of  current  at 
each  grid  point  were  factored  in  to  account  for  the  Coriolis  effect.  Measured  wind  data  from  a 
moored  buoy  in  the  vicinity  of  the  CODAR  coverage  area  were  used  to  calculate  the  covariance 
function  with  the  CODAR  data.  Winds  were  incorporated  in  two  ways.  First,  wind  observations 
prior  to  each  time  step  were  used  to  predict  future  surface  currents.  Next,  actual  wind 
observations  for  times  after  the  current  time  were  used  as  “predicted”  winds  to  evaluate  the 
extent  to  which  including  the  wind-driven  component  would  improve  performance. 

Of  these  enhancements,  the  addition  of  the  wind  contributed  to  a  marginal  improvement  in  the 
predictive  skill  of  the  STPS,  when  surface  current  fields  predicted  by  STPS  were  compared  to 
subsequent  CODAR  measurements.  The  incorporation  of  winds  into  STPS  produced  reductions 
in  RMS  prediction  error  ranging  from  2  cm/s  at  the  one  hour  forecast  time  to  less  than  1  cm/s  at 
prediction  times  of  a  few  hours,  or  approximately  10  percent.  The  weight  of  this  conclusion  is 
tempered  by  the  fact  that  winds  were  relatively  light  during  the  study  period.  The  impact  of  the 
wind  on  the  prediction  skill  may  be  underestimated  in  this  test,  and  further  study  may  be 
warranted  using  data  from  time  periods  when  winds  are  stronger  or  show  a  large  change  in 
direction  (e.g.  during  frontal  passages)  and  have  more  of  an  impact  on  surface  flow. 

7.2.  STPS  and  CODAR  prediction  of  SLDMB  trajectories 

Surface  currents  fields  derived  from  the  modified  STPS  and  the  CODAR  data  were  used  to 
generate  pseudo-trajectories  in  MAB  and  BIS.  The  small-scale  random  component  of  motion 
was  simulated  using  Monte-Carlo  RW  and  RF  models.  Trajectories  predicted  using  the  RW  and 
RF  models  were  compared  against  actual  SLDMB  trajectories.  The  method  used  to  calculate 
dispersion  was  found  to  affect  the  size  and  skill  of  the  RW  and  RF  approaches  relative  to  the 
SLDMB  trajectories.  The  magnitude  of  the  dispersion  coefficient  used  had  a  principal  effect  on 
the  search  area  size  and  the  success  of  the  RW  versus  RF  models. 

A  modification  not  envisioned  during  the  planning  of  this  project  was  to  incorporate  processing 
features  used  in  the  MAB  CODAR  processing  algorithm.  This  modification  to  STPS  included 
the  incorporation  of  Geometric  Information  System  (GDOP)  and  percentage  return  filtering, 
along  with  other  modifications.  The  incorporation  of  the  MAB  processing  algorithm  produced 
the  most  significant  improvement  in  the  skill  of  STPS,  reducing  the  RMS  error  of  the  STPS 
predictions  by  approximately  20  percent.  The  quality  of  the  returns  is  related  to  the  geometry  of 
the  transmitter  arrays  and  the  orientation  of  transmitters  relative  to  the  predominant  directions  of 
surface  current  flow.  These  results  may  indicate  that  future  efforts  to  improve  the  usability  of 
this  technology  should  also  examine  signal  processing  and  siting  issues.  In  summary,  the 
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greatest  improvement  to  the  usability  of  the  CODAR  system  resulted  from  an  improvement  to 
the  processing  of  CODAR  data  rather  than  from  adding  the  effect  of  forecast  winds  to  the  STPS. 

The  exercise  showed  that  CODAR  is  better  at  hindcasting  drifter  motion  in  these  areas  relative  to 
the  methods  currently  in  use.  The  most  significant  improvements  in  CODAR  skill  resulted  from 
CODAR  data  processing  and  filtering  steps.  Further  improvements  in  accuracy  may  be  realized 
in  improved  processing  of  returns  and  spatial  configuration  of  CODAR  networks.  The  time 
series  analysis  techniques  employed  in  this  study  produced  smaller  and  less  consistent 
improvements.  Modifications  to  the  time  series  analysis  approach  used  by  STPS  may  also  be 
beneficial. 

Future  improvements  to  the  predictive  performance  of  both  methods  will  also  require  more 
accurate  estimates  of  the  dispersion  coefficients  and  correlation  time  scale.  As  SLDMBs  are 
deployed  in  a  region,  current  velocity  data  could  be  saved  and  utilized,  in  conjunction  with 
CODAR  current  velocity  data,  to  improve  dispersion  parameter  estimates  in  the  region. 

7.3.  Factors  affecting  usability  of  HF  radar  products 

A  Monte-Carlo  simulation  task  was  undertaken  to  focus  on  improving  the  characterization  of 
uncertainties  in  the  CODAR  and  STPS  surface  current  fields.  Random  walk  (RW)  and  random 
flight  (RF)  Monte  Carlo  models  were  used  to  simulate  CODAR/STPS  error  and  dispersion  due  to 
turbulence  in  the  surface  current  field.  The  RW  and  RF  models  were  used  to  model  SLDMB 
trajectories  in  the  MAB  and  BIS.  Dispersion  coefficients  were  computed  using  variance  of 
CODAR  field  velocities  and  from  differences  between  CODAR  and  SLDMB  velocities. 
Dispersion  was  estimated  retrospectively  using  the  results  of  Monte  Carlo  simulations  of 
SLDMB  trajectories.  These  various  methods  resulted  in  estimates  of  the  dispersion  coefficient 
ranging  from  40  m2/s  to  700  m2/s. 

The  differences  between  search  areas  obtained  from  RW  and  RF  models  was  found  to  be  minor 
relative  to  the  differences  in  search  area  sizes  resulting  from  uncertainty  in  dispersion  (A).  The 
value  selected  for  the  dispersion  coefficient  (K)  emerges  as  the  critical  factor  in  defining  the  size 
and  growth  of  the  search  area  in  a  SAR  case,  and  subsequently  the  usefulness  of  the  CODAR  and 
STPS  products. 

7.4.  Operational  demonstration  of  SAROPS 

A  preliminary  version  of  SAROPS  EDS  was  used  to  collect  CODAR  current  data  and  STPS 
current  predictions  from  UCONN  and  Rutgers,  wind  observations  from  NDBC  buoy  sites,  and 
SLDMB  data  from  OSC  Martinsburg.  This  initial  test  of  the  retrieval  system  concluded  that 
remote  users  will  be  able  to  request  and  download  these  data  in  near-real  time  for  use  in  Search 
Area  Plans  (SAPs).  The  following  issues,  however,  need  to  be  considered  for  operational  use. 

•  The  SAR  operators  possessed  an  understanding  of  the  marine  environment  and  want 
more  information  on  the  data  that  they  use  for  SAR  cases.  They  need  information  on  the 
quality  and  reliability  of  the  environmental  data  that  they  use. 
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•  The  user  could  manually  select  data  sources  available  from  the  EDS,  but  the  system 
should  automatically  recommend  data  sources  and  parameters  such  as  regional 
dispersion/turbulence  based  on  the  AOI. 

A  30-day  evaluation  operational  test  and  evaluation  was  then  conducted  was  at  CG  Group 
Moriches  (GRUMOR).  SLDMBs  deployed  during  this  period  were  used  as  search  targets. 
Trajectory  predictions  were  made  using  STPS  and  CODAR.  Both  the  BIS  and  MAB  systems 
were  fully  functional  during  the  T&E  period.  Operators  expressed  interest  and  confidence  in  the 
CODAR/STPS  predictions.  When  SLDMB  trajectories  were  inside  the  CODAR/STPS  coverage 
area,  the  mean  value  and  variability  of  SLDMB-prediction  errors  using  CODAR/STPS  were 
measurably  lower  than  those  from  C2PC/JAWS.  When  the  drifter  was  outside  the 
CODAR/STPS  coverage  area,  the  mean  errors  were  comparable,  although  the  variability  of 
SLDMB-prediction  errors  obtained  from  CODAR/STPS  were  lower.  In  spite  of  this  result, 
caution  is  recommended  when  using  data  collected  outside  of  the  coverage  area. 

Significant  issues  that  emerged  were: 

•  Use  of  the  modified  SARMAP  to  interface  with  the  EDS  raised  issues  that  should  be 
considered  in  the  development  of  SAROPS,  e.g.,  decision  automation. 

•  The  lack  of  near  coastal  CODAR  coverage  was  a  concern  for  operators.  Future  versions 
of  the  EDS  will  need  to  blend  data  sources  such  as  CODAR  with  other  data  sources  to 
provide  complete  data  coverage  for  a  SAR  case. 

The  demonstration  project  had  no  method  to  assimilate  drift  data  from  the  SLDMB  tracks  into 
SARMAP.  The  tracks  could  be  visualized,  but  the  forecast  (and  hindcast)  solutions  did  not 
exploit  the  SLDMB  data.  A  system  could  be  established  to  catalog  SLDMB  and  CODAR  data 
for  individual  CODAR  regions.  Simply  maintaining  a  database  with  all  SLDMB  trajectories 
through  a  given  CODAR  region  along  with  concurrent  CODAR  velocities  will  allow  better 
estimates  of  prediction  error  and  dispersion  coefficients  for  a  specific  region.  The  necessary 
calculations  could  be  performed  on  a  regular  basis  (e.g.,  annually  or  seasonally). 

7.5.  Recommendations  for  the  U.S.  Coast  Guard 

CODAR  and  STPS  are  maturing  technologies  that  have  demonstrated  the  ability  to  improve  CG 
operational  planning.  Presently  however,  CODAR  and  STPS  are  not  ready  for  operational  use  in 
SAROPS.  The  following  recommendations  are  provided  to  make  them  viable  for  operational  use 
in  SAROPS. 

1)  Through  Oceans.US,  the  CG  should 

a.  Encourage  initiatives  to  improve  the  accuracy  and  coverage  of  HF  radar  arrays, 

b.  Encourage  development  of  enhancements  to  STPS  in  HF  radar  coverage  regions. 

c.  Partner  with  the  HF  radar  community  to  ensure  that  the  CG  requirements  are 
recognized. 


2)  In  the  SAROPS  developments  process,  G-OPR  should  include  methodologies  in  the  EDS  to 
fill  the  gaps  where  no  or  invalid  surface  current  data  exist.  This  applies  to  the  use  of  model 
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data  as  well  as  HF  radar  data.  In  areas  where  there  is  no  HF  radar  coverage,  alternative 
current  estimates,  such  as  historical  files  or  datum  marker  buoy  (DMB)  data,  should  be 
available. 

3)  In  conjunction  with  the  Operations  Systems  Center,  G-OPR  should  develop  an  archive  to 
store  SLDMB  data  and  periodically  calculate  and  store  regionally  representative  dispersion 
coefficients  from  SLDMB  trajectories.  Archived  SLDMB  data  could  also  be  matched  with 
seasonal  and  environmental  conditions  to  provide  better  historical  current  estimates  in  areas 
where  CODAR/STPS  is  not  available. 

4)  G-OPR  should  support  further  research  on  the  impact  of  forecast  wind  change  on  STPS 
current  estimates. 

7.6.  CODAR/HF  radar  community  specific  recommendations 

The  following  recommendations  are  provided  for  HF  radar  operators. 

1)  Add  standard-range  HF  radar  sites  to  fill  in  gaps  in  areas  closer  to  shore  and  in  approaches  to 
harbors,  especially  in  areas  where  spatial  variability  is  high.  Consider  site  development  by 
frequency  of  SAR  cases. 

2)  Document  the  software  that  is  used  to  process  data  in  each  area.  Meet  with  operators  in  other 
regions  to  develop  standards  for  combining  radials,  filling  in  missing  spatial  and  temporal 
data,  and  eliminating  low  quality  data.  Develop  data  sharing  standards  to  facilitate 
interfacing  with  STPS  and  SAROPS. 

3)  Develop  standard  quality  control  and  uncertainty  estimates  so  that  search  areas  can  be 
modeled  more  effectively. 

4)  Explore  schemes  to  use  stored  HF  radar  data  to  compute  and  regularly  update  regionally 
appropriate  dispersion  coefficients.  To  the  CG  user,  CODAR  surface  current  mapping 
credibility  is  contingent  on  reliable  updates  to  dispersion  coefficients. 
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APPENDIX  A. 

DESCRIPTION  OF  THE  RUTGERS  UNIVERSITY 
CODAR  DATA  PROCESSING 


A.l  Introduction 

High  frequency  radar  systems,  typically  deployed  along  the  coast  (Figure  A-l),  use  Bragg  peaks 
within  a  signal  (3  ~  30  MHz)  scattered  off  the  ocean  surface  to  calculate  radial  components  of 
the  total  surface  velocity  at  a  given  location  (Barrick,  Evens,  Weber,  1977).  Crombie  (1955) 
recognized  that  these  peaks  were  the  result  of  an  amplification  of  a  transmitted  wave  by  surface 
gravity  waves  with  a  wavelength  equal  to  half  that  of  the  transmitted  signal.  A  signal  scattered 
off  a  wave  and  back  toward  the  antenna  will  be  in  phase  with  a  signal  that  traveled  to  the  next 
surface  wave  (!4  transmit  wavelength  further)  and  returned  to  the  original  wave  (another 
Vi  transmit  wavelength).  The  frequency  of  the  backscattered  signal  will  be  shifted  depending  on 
the  velocity  of  the  scattering  surface.  Using  linear  wave  theory,  the  phase  speed  of  the  surface 
waves  can  be  separated  from  the  total  frequency  shift,  leaving  only  that  shift  due  to  the  surface 
current.  Over  a  given  time  period,  sites  along  the  coast  generate  radial  component  maps  of  the 
surface  current  with  typical  resolutions  on  the  order  of  1  to  6  km  in  range  and  5  degrees  in 
azimuth  (Figure  A-2).  Because  the  Doppler  shift  can  only  resolve  the  component  of  the  current 
moving  toward  or  away  from  the  site,  information  from  at  least  two  sites  must  be  geometrically 
combined  to  generate  total  surface  current  maps. 


Figure  A-l.  Transmit  (left)  and  receive  (right)  antennas  for  a  typical  CODAR-type  long-range 
system. 
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Figure  A-2.  Sample  radial  vector  coverage  for  a  long-range  and  standard-range  HF  radar  system. 


Rutgers  University  operates  an  array  of  10  CODAR-type  HF  radar  systems,  6  of  which  are  lower 
resolution,  long-range  sites.  With  an  operating  frequency  around  5  MHz,  these  sites  measure 
surface  currents  within  the  upper  2.4  meters  of  the  water  column  (Stewart  and  Joy,  1974). 
Typical  spatial  resolutions  are  on  the  order  of  6  km,  with  maximum  ranges  exceeding  200  km. 
Four  sites  along  the  coast  of  New  Jersey  from  Wildwood  to  Sandy  Hook  provide  hourly  surface 
current  maps  over  the  New  Jersey  Shelf  (Figure  A-3).  These  four  sites  form  one  cluster  of 
systems  within  the  North  East  Observing  System  (NEOS).  These  sites  use  global  positioning 
system  (GPS)  synchronization  so  that  each  site  is  operating  at  the  same  frequency  and 
bistatically  linked  to  the  other  sites  within  the  network. 


Figure  A-3.  Surface  current  map  of  January  19,  2005  at  20:00  GMT.  Coverage  is  shown  for  the 
long-range  (left)  and  standard-range  (right)  systems. 
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Two  standard-range  systems,  currently  deployed  on  opposite  sides  of  the  entrance  to  New  York 
Harbor  (Figure  A-3),  are  nested  within  the  New  Jersey  long-range  cluster.  With  an  operating 
frequency  of  25  MHz,  these  systems  measure  the  current  within  the  upper  50  cm  of  the  water 
column.  Typical  spatial  resolutions  are  on  the  order  of  1  km  with  maximum  ranges  out  to  40  km. 
The  two  sites  were  recently  moved  from  their  original  location  along  the  southern  coast  of  New 
Jersey  to  the  New  York  Bight  apex  to  support  (1)  a  National  Science  Foundation  (NSF)-funded 
project  focusing  on  the  Hudson  River  outflow  and  (2)  the  development  of  ship-tracking 
capabilities  near  the  harbor  entrance.  On  26  January  2005,  a  third  standard  range  site  operated 
by  Stevens  Institute  was  deployed  within  the  lower  harbor.  This  site  paired  with  the  two  Rutgers 
sites  brings  existing  high-resolution  coverage  into  the  harbor. 

A.2  Radial  data  processing 

These  CODAR-type  systems  are  direction-finding  systems  that  use  a  three-element  receive 
antenna  mounted  on  a  single  post  to  determine  the  direction  of  the  incoming  signals.  The 
angular  resolution,  set  in  the  processing,  is  typically  5  degrees  (Teague,  Vesecky,  Fernandez, 
1997;  Barrick  and  Lipa,  1996).  For  the  radial  vector  map  processing,  raw  time  series  from  each 
element  are  transformed  into  cross-spectra  with  a  1024-point  fast  Fourier  transform  (FFT).  For 
the  long-range  systems,  these  data  are  created  every  17  minutes  and  averaged  hourly,  centered  on 
each  half  hour.  For  each  range  bin,  spectra  are  run  through  a  multiple  signal  classification 
(MUSIC)  algorithm  with  the  measured  beam  patterns  (Kohut  and  Glenn,  2003)  to  calculate  the 
bearing  of  each  radial  velocity  within  the  first-order  Bragg  region  (Schmidt,  1986;  Barrick  and 
Lipa,  1999).  This  first-order  region  was  limited  to  signal  to  noise  ratios  (SNR)  of  at  least  5.0  dB. 
Each  hour,  seven  of  these  hourly  velocity  estimates,  each  overlapped  by  half  an  hour,  are 
averaged  to  get  a  4-hour  running  average  for  the  center  time.  This  file  is  input  to  the  vector 
combination  procedure  to  generate  the  total  vector  file  described  in  section  A.3. 

The  number  of  radial  vectors  in  any  given  range  cell  depends  on  (1)  the  digitization  interval 
determined  by  the  FFT  length,  operating  frequency,  and  sweep  rate  and  (2)  the  number  of 
antenna  elements  in  the  receive  array.  The  digitization  interval  controls  how  many  radial 
velocity  vectors  are  available,  and  the  number  of  antenna  elements  determines  the  number  of 
possible  MUSIC  solutions  for  each  radial  velocity.  For  our  setup  with  three-receive  antenna 
elements,  a  radial  velocity  every  3.22  cm/s  can  be  placed  in  up  to  two  angular  bins. 

Consequently,  if  there  are  more  than  two  angles  with  a  given  radial  velocity  or  periods  of  weak 
surface  currents,  data  coverage  will  be  reduced.  Based  on  these  constraints,  an  interpolation  was 
done  in  azimuth.  The  scheme  linearly  interpolated  radial  data  across  angular  gaps  of  15  degrees 
or  less  (Kohut,  Roarty,  Glenn,  2004;  Kohut  and  Glenn,  2003). 

A.3  Total  vector  combination 

Every  hour,  the  available  radial  vector  maps  were  geometrically  combined  into  a  single  total 
vector  map.  All  radial  component  vectors  within  10  km  of  each  grid  point  were  used  in  the 
combination.  A  total  vector  was  generated  only  if  at  least  three  radial  vectors  from  at  least  two 
remote  sites  were  used  in  the  combination.  When  the  orthogonal  velocity  components  are 
combined  into  total  vector  estimates,  errors  are  introduced.  The  uncertainty  of  the  combined 
totals  can  be  separated  into  radial  vector  uncertainty  and  geometric  uncertainty.  The  geometric 
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uncertainty  is  based  on  the  angles  of  the  radial  component  vectors.  The  farther  the  radials  are 
from  orthogonal,  the  greater  the  uncertainty.  This  uncertainty  is  described  by  Chapman  and 
Graber  (1997)  as  geometric  dilution  of  precision  (GDOP).  These  uncertainties  increase  along  the 
baseline  between  the  coastal  sites  as  well  as  offshore.  In  these  regions,  the  radial  component 
vectors  are  resolving  the  along-shore  and  cross-shore  velocity  components,  respectively.  The 
relative  magnitude  of  the  geometric  uncertainty  was  calculated  at  each  grid  point  for  every 
current  map  (Gurgel,  1994).  Using  this  scalar  as  an  indicator  of  the  magnitude  of  the  geometric 
contribution  to  the  uncertainty,  data  subject  to  poor  geometry  were  eliminated.  This  was  done  by 
eliminating  all  data  above  a  threshold  of  1.25.  This  1.25  threshold  was  chosen  based  on 
qualitative  analysis  of  previous  data.  Drifter  deployments  offer  an  excellent  opportunity  to  better 
define  this  threshold  and  apply  it  to  future  system  operation.  For  the  duration  of  the  project,  total 
vector  fields  were  made  with  and  without  the  geometric  filter  (Figure  A-4).  For  the  purpose  of 
this  analysis,  the  unfiltered  data  was  used  to  evaluate  the  filter,  and  the  filtered  data  was  used  for 
all  direct  comparisons  to  drifter  velocities. 
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Figure  A-4.  Surface  current  map  combined  with  (right)  and  without  (left)  the  geometric  filter. 


A.4  Drifter  velocities 

The  drifter  position  data  were  used  to  evaluate  the  CODAR  observations.  The  drifter  data  used 
here  are  from  the  deployment  during  the  summer  of  2004  (July  27th  to  August  3 1st).  Throughout 
the  deployment,  drifter  positions  were  reported  every  half  hour  (Figure  A-5).  SLDMBs  #32773, 
#32776,  and  #32779  were  deployed  in  the  BIS  and  moved  southwest  along  the  southern  Long 
Island  coast.  SLDMBs  #43057  and  #43060  were  deployed  in  the  upper  Hudson  Canyon  and 
moved  south  into  the  center  of  the  coverage,  and  SLDMBs  #43061  and  #43062  were  deployed 
farther  offshore  and  moved  in  and  out  of  the  coverage  through  the  deployment.  Drifter  velocities 
were  calculated  to  match  the  averaging  of  the  CODAR  processing.  Velocities  based  on  two 
drifter  positions  one  hour  apart  were  calculated  every  half  hour.  Each  hour  seven  of  these  hourly 
velocity  estimates,  each  overlapped  by  half  an  hour,  were  averaged  to  get  a  4-hour  running 
average  for  the  center  time.  These  drifter  data  were  compared  to  CODAR  observations.  At  each 
time  step,  the  CODAR  data  were  linearly  interpolated  to  the  drifter  position.  Because  only  radial 
data  within  10  km  of  a  grid  point  were  combined  to  create  the  total  vector  grid,  comparisons 
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were  limited  to  times  when  the  drifter  was  no  more  than  10  km  from  the  closest  CODAR  grid 
point. 


Longitude  (DegreesMnutes) 
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Figure  A-5.  Tracks  of  the  seven  drifters  (black  line)  deployed  between  July  27,  2004  and 
August  31,  2004.  The  deployment  sites  (green  circle)  and  the  50  percent 
CODAR  coverage  contour  (gray  line)  are  also  shown. 

A.5  CODAR/drifter  comparisons 

CODAR  data  were  validated  against  seven  drifters  deployed  in  the  summer  of  2004.  For  each 
drifter  within  the  CODAR  coverage  (Figure  A-5),  comparisons  were  drawn  between  spatially 
interpolated  CODAR  observations  and  climatological  estimates.  Because  the  CODAR  data 
sampled  the  surface  current  over  different  horizontal  and  spatial  scales  than  the  drifter,  the 
current  estimates  differed  in  the  absence  of  instrument  uncertainty.  Therefore,  the  contribution 
to  the  RMS  difference  includes  CODAR  uncertainty,  drifter  uncertainty,  and  sub-grid  scale 
uncertainty.  The  standard  deviation  of  the  drifter  data  calculated  with  the  seven  half-hourly  data 
averaged  each  hour  was  used  as  an  estimate  of  the  sub-grid  scale  uncertainty.  For  reference,  the 
average  standard  deviation  for  all  drifters  through  the  entire  deployment  was  7.2  cm/s  for  the  east 
component  and  6.2  cm/s  for  the  north  component. 

The  average  RMS  difference  among  all  available  drifter  and  CODAR  velocities  was  9. 1  cm/s  in 
the  east  component  and  12.1  cm/s  in  the  north  component.  The  same  drifter  data  were  compared 
to  Mariano  climatology  giving  an  RMS  difference  of  17.3  cm/s  in  the  east  direction  and 
15.9  cm/s  in  the  north  direction  (Table  A-l). 
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Table  A-l.  Comparison  statistics  among  drifter  velocity,  Mariano  climatology,  and  CODAR 
observations.  For  reference,  the  standard  deviation  of  the  drifter  data  was  7.2  cm/s 


for  the  easi 

direction  and  6.2  cm/s  for  the  north  component. 

SLDMB 

Number 

Velocity 

Component 

SLDMB/CODAR 
RMS  (cm/s) 

SLDMB/ 

Mariano  Climatology 
RMS  (cm/s) 

32773 

East 

7.7 

18.1 

32773 

North 

11.0 

13.5 

32776 

East 

10.1 

19.4 

32776 

North 

12.4 

15.3 

32779 

East 

7.6 

17.3 

32779 

North 

10.4 

13.3 

43057 

East 

8.1 

16.9 

43057 

North 

11.3 

14.9 

43060 

East 

9.6 

16.6 

43060 

North 

10.1 

15.2 

43061 

East 

12.2 

17.9 

43061 

North 

16.9 

22.7 

43062 

East 

8.6 

15.2 

43062 

North 

12.6 

16.1 

Average 

East 

9.1 

17.3 

Average 

North 

12.1 

15.9 

For  each  drifter,  the  RMS  difference  with  the  CODAR  observations  was  lower  than  that  with 
Mariano  climatology.  In  addition,  the  east  component  of  the  CODAR  data  had  a  lower  RMS 
difference  than  the  north  component  for  all  drifters.  This  difference  is  evidence  of  the  additional 
uncertainty  in  the  north  component  due  to  site  geometry.  All  the  sites  in  this  network  were 
deployed  along  the  New  Jersey  coast.  Given  the  drifter  location  and  the  radial  coverage  of  these 
sites,  the  east  component  of  the  total  vectors  was  better  resolved  by  the  radial  vectors  than  the 
north  component.  SLDMB  #32773  clearly  illustrated  the  geometric  uncertainty  in  the  vector 
combination  (Figure  A-6).  Early  on  in  the  deployment,  the  drifter  was  south  of  Long  Island 
where  the  radial  vectors  measured  at  the  CODAR  sites  better  resolved  the  east  velocity 
component.  The  time  series  clearly  shows  better  agreement  between  the  east  component  of  the 
CODAR  and  drifter  observations  compared  to  the  north  component.  As  the  drifter  moved  closer 
to  the  sites  and  better  geometry,  the  north  component  was  better  resolved  and  the  agreement 
improved.  SLDMB  #43060  spent  more  time  in  the  coverage  area  with  better  geometry,  and  the 
agreement  was  good  throughout  the  time  series.  The  CODAR  data  are  an  improvement  to 
existing  Mariano  climatology  when  compared  to  drifter  velocity  estimates.  Future  work  will 
focus  on  refining  the  geometric  filter  and  processing  algorithms  to  improve  both  the  geometric 
and  CODAR  uncertainties. 
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Figure  A-6.  Drifter  track  (red  line)  (left),  CODAR  coverage  area  (green)  (left),  and  time  series 
of  the  east  and  north  velocity  components  from  drifter  locations  (red  trace), 
CODAR  (blue  trace),  and  Mariano  climatology  (black  line).  The  upper  panel  is  for 
SLDMB  #32773,  and  the  lower  panel  is  for  SLDMB  #43060. 
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APPENDIX  B. 

STPS  ENHANCEMENTS  IN  THE  BIS  REGION 


B.l  Enhanced  STPS  with  no  wind  in  the  BIS  region 

Using  CODAR  data  for  the  December  2003  to  January  2004  period,  the  period  when  drifter  data 
were  available  in  the  BIS  area,  there  was  no  substantial  improvement  to  the  performance  of  the 
algorithm  when  local  covariance  estimates  and  the  covariance  of  the  current  components  were 
used.  Figure  B-l  shows  a  comparison  of  the  area-averaged  error  in  the  forecasts  of  the  original 
algorithm  with  that  of  the  new  algorithm  in  support  of  this  conclusion.  The  difference  in  root- 
mean-squared  (RMS)  errors  is  less  than  0.25  cm/s. 


Figure  B-l.  Area-averaged  RMS  error  in  the  east  (blue  lines)  and  north  (red  lines)  components 
of  the  surface  currents  estimated  by  the  original  (solid  lines)  STPS  algorithm  and 
the  enhanced  STPS  algorithm  with  local  covariance  and  covariance  component 
modifications  (dashed  lines)  in  the  BIS  region. 


The  spatial  structure  of  the  errors  at  the  6-hour  forecast  time  is  presented  in  Figure  B-2.  The 
distribution  and  magnitude  of  the  errors  are  very  similar  in  large  areas  where  the  error  is 
approximately  10  cm/s  and  areas  of  slightly  higher  error  where  the  geometry  of  the  CODAR 
system  degrades  the  precision  of  the  vector  estimates.  The  top  figures,  the  east  components, 
have  more  areas  of  green  and  yellow,  corresponding  to  higher  RMS  errors. 
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Figure  B-2.  The  distribution  of  the  time-averaged  RMS  error  in  the  surface  current  predictions 
of  original  (left  two  frames)  and  new  STPS  (right  frames)  in  the  BIS  area.  The 
errors  in  the  east  components  are  presented  in  the  top  two  frames  and  the  north 
component  in  the  lower  frames. 


B.2  Enhanced  STPS  with  measured  wind  in  the  BIS  region 

Measured  wind  data  were  included  in  the  Gauss-Markov  procedure  to  generate  current  forecasts. 
Note  that  wind  forecasts  were  not  included  at  this  stage  of  the  analysis.  The  wind  data  were 
acquired  from  NDBC  buoy  #44025.  Figure  B-3  displays  a  comparison  of  the  area-averaged 
error.  Again,  there  was  no  improvement  in  the  prediction  skill. 
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Figure  B-3.  Area-averaged  RMS  error  in  the  predicted  east  (blue  lines)  and  north  (red  lines) 
velocity  components  of  the  surface  currents  estimated  by  the  enhanced  STPS 
algorithm  without  wind  data  (dashed  lines)  and  the  enhanced  algorithm  with 
measured  wind  data  (solid  lines)  in  the  BIS  area. 
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Abstract 


The  ability  to  predict  the  movement  of  particles  drifting  at  the  surface  of  the  coastal  ocean  is  of 
great  value  to  agencies  responsible  for  the  management  of  oil  spills,  search  and  rescue 
operations,  marine  safety  and  coastal  water  quality.  High  frequency  radio  surface-current 
measuring  systems  that  are  capable  of  mapping  wide  areas  can  be  exploited  for  this  purpose.  We 
present  a  simple  statistical  technique  that  exploits  observations  to  extrapolate  Eulerian  velocity 
estimate  12  hours  into  the  future.  We  then  use  these  predictions  in  the  Euler-Lagrange 
transformation  to  simulate  the  trajectories  of  drifting  particle.  The  performance  of  both  the 
Eulerian  forecast  and  the  trajectory  forecast  technique  are  evaluated  by  comparison  to 
observations  in  Block  Island  Sound  and  the  adjacent  shelf.  We  also  evaluate  the  performance  to 
an  algorithm  similar  to  that  in  operational  use  to  gauge  the  practical  value  of  the  system.  We 
find  that  the  surface  current  predictions  have  root  mean  square  errors  averaged  over  the  -800 
km2  area  of  less  than  15  cm/s  for  forecasts  up  to  12  hours.  The  errors  appear  to  be  correlated 
with  wind,  suggesting  that  the  algorithm  could  be  further  refined.  When  employed  to  forecast 
drifter  trajectories  we  find  that  the  root  mean  squared  separation  between  predicted  drifter 
locations  and  observations  show  that  95  percent  of  the  differences  are  less  that  7  km  at  6  hours 
and  13  km  at  12  hours.  This  is  approximately  a  factor  of  two  better  than  methods  currently  in 
operational  use. 

C.l  Introduction 

Several  surface  current  mapping  systems  (e.g.  CODAR,  WERA,  OSCAR)  that  exploit  the 
Doppler  shift  in  the  Bragg  scatter  from  high  frequency  radio  transmissions  have  been  developed 
and  are  now  commercially  available.  We  will  refer  to  them  collectively  as  SCMS.  Though  they 
differ  in  important  technical  characteristics,  we  believe  our  work  is  equally  applicable  to 
observations  from  all  designs.  Some  installations  have  been  operating  for  more  than  7  years 
(Glenn  and  Schofield,  2003)  and  the  data  has  been  extensively  used  in  science  programs  (Kohut 
et  al,  2004).  Though  there  is  a  variety  of  algorithms  that  transform  the  backscatter  spectra  to 
velocity  and  many  technical  improvements  remain  possible,  system  designs  have  become  stable 
and  very  reliable  in  the  last  few  years  and  are  ready  for  operational  use  as  a  central  component  of 
an  integrated  ocean  observing  system. 
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The  timely  and  accurate  prediction  of  the  motion  of  material  floating  at  the  ocean  surface  is  a 
critical  limitation  to  the  effective  response  to  marine  safety  emergencies  and  the  accidental 
release  of  pollutants.  In  particular,  the  availability  of  SCMS  has  the  potential  to  offer  substantial 
advantage  to  the  CG,  which  currently  relies  on  sparse  network  stations  at  which  tidal  current 
harmonics  have  been  estimated.  As  has  been  well  established  (see,  for  example,  Ridderinkof  and 
Zimmerman,  1992),  the  interaction  of  tidal  currents  with  spatially  varying  mean  flows  can  create 
extremely  rapid  Lagrangian  dispersion.  Since  SCMS  can  reveal  the  spatial  structure  of  the 
circulation,  they  offer  the  potential  of  substantial  improvement  in  the  effectiveness  of  search 
planning.  Because  of  the  shape  of  the  coastline,  harbors  and  estuaries  have  complex  currents,  a 
lot  of  boat  traffic,  and  frequent  marine  accidents.  SCMS  could,  therefore,  help  save  lives  in  the 
immediate  future. 

The  spectral  averaging,  data  telemetry  and  processing  required  by  SCMS  delays  the  availability 
of  observations.  The  most  recent  observations  currently  available  can  be  1  to  3  hours  old.  The 
integration  of  surface  current  measuring  infrastructure  to  operational  safety  and  management 
systems,  therefore,  requires  that  a  forecast  capability  be  developed  for  the  Eulerian  currents  and 
that  an  algorithm  be  developed  that  can  effectively  exploit  forecasts  to  predict  Lagrangian 
trajectories  (particle  paths). 

In  this  paper  our  intent  is  to  describe  the  development  and  performance  of  an  algorithm  that 
exploits  data  to  produce  short  term  (0-24  hour)  forecasts.  The  approach  we  adopted  makes  no 
assumptions  about  the  underlying  dynamics  and  ignores  the  coastal  geometry  and  bathymetry. 
The  primary  advantage  of  this  approach  is  that  the  algorithm  can  be  readily  applied  in  any  area 
with  an  operational  SCMS  with  limited  effort  and  expense.  It  is  unlikely  that  this  will  be  the  best 
approach  in  the  long  run.  However,  it  is  easily  implemented  and  it  establishes  a  benchmark 
against  which  more  sophisticated  algorithms  can  be  evaluated. 

We  use  the  Eulerian  forecast  algorithm  to  predict  Lagrangian  trajectories  in  Block  Island  Sound 
and  the  adjacent  continental  shelf.  To  evaluate  the  effectiveness  of  the  approach  we  deployed 
surface  drifters.  To  establish  the  relative  importance  of  the  limitation  imposed  by  the  data  and 
the  prediction  system  we  simulated  the  drifter  trajectories  using  data  and  the  forecast  system. 

We  also  evaluated  the  performance  of  the  system  that  is  currently  in  use  that  neglects  much  of 
the  spatial  structure  in  the  flow  field.  We  demonstrate  that  the  SCMS-based  system  is  limited  by 
data  quality  but  is  substantially  better  than  the  available  technology  in  the  test  area. 

In  the  next  section  we  describe  the  characteristics  of  the  flow  in  the  area  of  Block  Island  Sound 
where  a  standard  range  CODAR  system  has  been  in  operation  since  2000.  We  deployed  a 
flotilla  of  GPS  tracked  CODE-style  drifters  to  assess  the  effectiveness  of  the  system  and  describe 
this  data  in  Section  3.  We  then  define  the  empirical  forecasting  algorithm  in  Section  4  and  assess 
its  performance  by  comparing  predictions  to  SCMS  observations  in  Section  5.  The  effectiveness 
of  the  simulation  of  drifter  trajectories  using  both  the  observed  and  forecast  currents  is  assessed 
in  Section  6  by  comparing  them  to  drifter  tracks.  We  summarize  the  results  and  comment  on 
potential  developments  in  Section  7. 
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C.2  Block  Island  Sound  current  observations 


In  2000  the  National  Ocean  Partnership  Program  (NOPP)  sponsored  the  development  of  an 
experimental  coastal  ocean  observatory  at  the  approaches  to  Block  Island  Sound,  an  area  of 
strong  tides  and  significant  buoyancy  driven  circulation  on  the  Southern  New  England  Shelf. 
Figure  C-l.  The  coastline  and  bathymetry  of  Block  Island  Sound  with  the  location  and 
approximate  observation  areas  of  three  CODAR  sites,  shows  a  map  of  the  area  with  the 
bathymetry.  Ullman  and  Comillion’s  (1999)  analysis  of  surface  temperature  observations 
demonstrated  that  this  is  an  area  with  persistent  and  complicated  fronts.  O’Donnell  and  Allen 
(1992)  used  drifter  observations  and  data  from  a  moored  current  meter  to  demonstrate  that  there 
is  substantial  mean  exchange  flow  through  Block  Island  Sound  associated  with  the  estuarine 
outflow  from  Long  Island  Sound.  They  showed  the  surface  motion  is  to  the  east  and  south  at 
approximately  10  cm/s  carrying  the  brackish  effluent  to  the  adjacent  continental  shelf.  Salty 
shelf  water  flows  north  and  west  into  Long  Island  Sound  along  the  relic  river  channel  that  is 
apparent  in  the  isobaths. 
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Figure  C-l.  The  coastline  and  bathymetry  of  Block  Island  Sound  with  the  location  and 
approximate  observation  areas  of  three  CODAR  sites. 

Figure  C-l  also  shows  the  location  of  the  three  standard  range  CODAR  SeaSondes®  that  were 
deployed  in  2000  at  Misquamicut,  RI,  Montauk  Point,  NY,  and  Block  Island,  RI,  also  displayed 
with  their  approximate  areas  of  coverage.  HF  radio  frequency  backscatter  observations  by  these 
instruments  can  be  combined  to  obtain  estimates  of  the  near  surface,  approximately  the  top  2 
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meters,  velocity  component  in  the  radial  direction  towards  the  receivers  in  2  km  range  and  5 
degree  azimuth  cells  every  hour.  At  the  FRONT  site  the  “radials”  from  all  three  sites  are 
combined  by  weighted  area  averaging  to  form  east  and  north  components  of  the  velocity  at  the 
points  indicated  by  the  “+”  symbols  in  Figure  C-l.  Ullman  and  Codiga  (2004)  present  an 
analysis  of  the  mean  and  seasonal  variation  of  the  circulation  in  the  area  using  observations  from 
this  system. 

C.3  Drifter  Data 

The  RV  Connecticut  deployed  four  CODE-style  surface  drifters  (see  Davis,  1985)  equipped  with 
GPS  navigation  and  Argos  data  telemetry  systems  in  Block  Island  Sound  during  December  2002 
within  the  coverage  area  of  the  FRONT  CODAR  system.  The  four  rectangular  drag  elements 
(sails)  of  the  drifter  were  0.7m  by  0.3m  and  were  suspended  0.3m  below  the  surface.  Several  of 
these  drifters  were  recovered  and  redeployed,  while  one  failed  to  return  useful  data. 
Subsequently,  in  March  2003,  CG  helicopters  air-deployed  an  additional  12  drifters,  of  which  4 
were  deployed  in  the  FRONT  region  and  8  off  the  New  Jersey  shore  for  a  study  that  will  be 
described  elsewhere.  The  drifters  reported  their  positions  at  0.5-hour  intervals  via  the  Argos 
network.  The  time  series  of  drifter  position  were  decimated  to  hourly  intervals  to  align  with  the 
hourly  CODAR  current  estimates.  Each  drifter  track  was  broken  up  into  a  series  of  24-hour 
segments  that  were  treated  as  independent  trajectories  for  the  analysis  of  errors  in  drift  prediction 
performance  (Figure  C-2).  The  start  times  of  the  segments  were  offset  by  12  hours.  Thus  the 
first  trajectory  for  a  given  drifter  starts  at  the  time  of  the  first  good  position  (to)  with  subsequent 
pseudo-independent  trajectories  starting  at  to+ 12,  to+ 24,  etc. 
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Figure  C-2.  Trajectories  of  drifters  launched  in  the  area  of  Block  Island  Sound  in  December 
2002  and  March  2003. 

C.4  Eulerian  Current  Forecast  Algorithm 

C.4.1  Overview 

The  surface  current  forecast  algorithm  is  based  on  an  empirical  decomposition  of  the  currents 
into  three  components.  As  in  most  of  the  coastal  ocean,  a  strong  component  of  the  motion  is  due 
to  tides,  a  weaker  and  more  slowly  varying  component  is  driven  by  winds,  and  a  more  persistent 
motion  associated  with  ocean  density  variations  and  very  low  frequency  waves.  The  tidal 
circulation  is  periodic  and,  therefore,  easy  to  forecast.  The  largest  five  tidal  constituents  in  the 
FRONT  region  are  listed  in  Table  C-l.  Three  are  semidiurnal  and  two  are  diurnal.  Least- 
squares  harmonic  analysis  method  is  used  to  compute  the  amplitude  and  phase  of  these 
constituents  at  each  CODAR  grid  point.  In  our  forecast  algorithm  the  harmonic  analysis  is 
performed  each  week  using  the  observations  acquired  during  the  prior  month.  The  constituents 
are  used  to  predict  the  tidal  component  of  the  circulation  for  the  next  week. 

C.4.2  Tidal  Circulation 

The  least-squares  method  for  determination  of  the  harmonic  constituent  amplitudes  and  phases 
also  provides  estimated  uncertainties  in  these  parameters  (Press,  et  al.,  1992).  Using  standard 
propagation  of  errors  techniques  (Emery  and  Thomson,  1997)  these  uncertainties  can  be  used  to 
estimate  the  error  in  a  tidal  prediction  made  with  the  harmonic  constituents.  The  prediction  errors 
are  independent  of  time,  but  exhibit  some  variation  with  location  as  shown  in  Figure  C-3. 
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Generally,  the  errors  are  less  than  3-4  cm/s.  In  a  small  region  southeast  of  Block  Island,  however, 
they  reach  5-10  cm/s. 


Table  C-l.  Major  tidal  constituents  in  FRONT  area. 


Constituent 

Period  (hours) 

s2 

12.00 

m2 

12.42 

n2 

12.66 

Ki 

23.93 

0! 

25.82 
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U  (East)  Tidal  Prediction  Uncertainty 


V  (North)  Tidal  Prediction  Uncertainty 
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Figure  C-3.  Uncertainty  in  the  eastward  (top)  and  northward  (bottom)  components  of  the  tidal 
current  predictions.  The  predictions  are  made  using  constituents  derived  from 
harmonic  analysis  of  1 -month  time  series. 
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C.4.3  Non-tidal  Motion:  Hedging 


Two  methods  for  predicting  the  non-tidal  current  component  were  tested.  We  refer  to  the  first  as 
the  “hedging”  method.  In  it  we  assume  that  the  non-tidal  current  observed  for  the  last  TH  persists 
throughout  the  prediction  interval.  The  second  method  exploits  the  archived  data  to  characterize 
the  structure  of  the  correlation  in  the  data.  It  is  essentially  an  application  of  Gauss-Markov 
estimation  (Wunsch,  1996). 


The  hedging  method  requires  the  estimation  of  the  non-tidal  current,  u'NT ,  at  the  time  of 
prediction,  t  .  The  simplest  approach  is  to  subtract  the  predicted  tidal  current,  uT ,  from  the 
observed  current,  uOBS ,  at  the  time  of  the  prediction.  However,  to  reduce  noise,  some  averaging 
is  desirable.  Our  hedging  estimate  of  the  non-tidal  flow  is  computed  as 

UNt(J  p)  =  NT^f  p-m)  (1) 

^  m= 1 

where  u'NT  =  uOBS  -  uT  and  M  is  the  number  of  samples  in  the  averaging  interval  TH  .  The 
choice  of  TH  is  somewhat  arbitrary  but  should  reflect  the  timescale  of  autocorrelation.  We  have 
tested  the  values  of  25  and  49  hours.  The  predictions  derived  from  the  25-hour  average  were 
more  accurate  for  short-term  predictions  (i.e.  24  hours  or  less)  and  we  used  it  in  the  operational 
evaluation  discussed  below. 

C.4.4  Non-tidal  Motion:  Gauss-Markov 

It  has  been  well  established  that  the  non-tidal  motion  in  this  area  is  autocorrelated  (Pettigrew, 
1981)  as  it  is  elsewhere  in  the  coastal  ocean  (see,  for  example,  Kundu  and  Allen,  1976).  This 
correlation  can  be  exploited  by  representing  the  desired  prediction  at  time  tp+ St  as  a  linear 
combination  of  previous  measurements,  i.e. 

M 

UNT  (?p  +  =  ^  amUNT  ft  p-m  )  >  (2) 

m=\ 


where  am  represent  weights  and  the  sum  is  taken  over  the  M  prior  data  values.  Unlike  the 
hedging  method  defined  in  equation  (1),  the  predicted  non-tidal  current  is  not  constant  during  the 
prediction  interval  since  (2)  is  evaluated  for  each  tp+St.  The  Gauss-Markov  theorem  states  that 
appropriately  chosen  weights,  based  on  knowledge  of  both  the  measurement  errors  and  the 
autocovariance  of  the  data,  produce  a  solution  that  minimizes  the  variance  of  the  difference 
between  the  estimate  and  the  true  value.  For  a  data  series  with  zero-mean,  this  estimate  is 
unbiased  and  is  referred  to  as  the  best  linear  unbiased  estimate  (Wunsch,  1996). 

Our  procedure  assumes  that  measurement  errors  have  zero  mean,  are  uncorrelated,  and  are 
normally  distributed  with  standard  deviation  distributed  as  in  the  maps  of  Figure  C-4. 
Autocovariance  functions,  are  calculated  at  each  location  from  one  month  of  data  and  averaged 
over  the  entire  CODAR  domain.  The  autocovariance  functions  ( Rxx ),  shown  in  Figure  C-5, 
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describe  the  degree  of  temporal  correlation  in  the  data  time  series, 
weights  are  found  to  be 

«/ = z  k*  (tP ,  tj  )K  it,  ,tj)+ Rnn  it, ,  tj  )]_1 }. 

j= i 


If  we  define  Rnn=<jl,  the 

(3) 
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Figure  C-4.  Estimated  uncertainty  in  CODAR  current  observations. 
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U  -  normalized  autocovariance  function 


Figure  C-5.  Autocovariance  functions  for  u  (top)  and  v  (bottom),  averaged  over  the  entire 
CODAR  domain  for  the  month  of  November  2002. 

Data  statistics  are  assumed  to  be  stationary  in  time,  but  in  operational  mode  this  calculation 
could  be  repeated  regularly,  perhaps  once  a  month,  to  ensure  that  long-term  variations  of  the 
statistical  character  of  the  fields  are  captured. 

An  attractive  advantage  of  the  Gauss-Markov  approach  is  that  it  also  produces  estimates  of  the 
uncertainty  in  the  prediction.  At  short  time  scales,  and  in  cases  where  the  measurement  signal  to 
noise  ratio  is  high,  this  uncertainty  is  approximately  equal  to  the  measurement  noise.  As  the  time 
of  the  prediction  extends,  the  uncertainty  approaches  the  data  variance.  This  estimate  of  the 
uncertainty  can  be  very  useful  for  a  short-term  predictive  modeling  system,  and  its  consistency  is 
evaluated  below. 

C.5  Evaluation  of  Forecasts 

The  algorithm  is  evaluated  by  comparing  the  predicted  surface  currents  to  subsequent 
measurements.  The  two  methodologies  are  tested  using  a  28-day  period  starting  December  1, 
2002.  The  measure  of  the  prediction  error  is  the  root-mean-square  (RMS)  difference  between 
predicted  and  observed  currents,  defined  for  the  eastward  component  as: 


where  up(tn)  is  the  predicted  eastward  current,  and  n0bs(tn)  is  the  corresponding  measurement  at 
time  tn,  and  the  average  is  spatial  over  N  pairs  of  observations/predictions.  URMS  is  expected  to 

depend  on  the  time  since  the  prediction  was  made,  tn  - t  ,  the  forecast  lag.  Equation  (4)  and  the 
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corresponding  equation  for  the  northward  velocity  component,  are  evaluated  separately  for  each 
forecast  lag  up  to  a  maximum  lag  of  25  hours. 

C.5.1  Region-Averaged  Differences 

Evaluation  of  the  differences  averaged  over  the  entire  FRONT  domain  provides  an  assessment  of 
the  performance  of  the  algorithm.  Figure  C-6  shows  the  dependence  of  Ums  for  both  the 

hedging  and  Gauss-Markov  methods  functions  of  the  forecast  lag.  Note  that  the  same  tidal 
prediction  method  is  used  for  both  cases.  As  expected,  the  errors  increase  with  the  forecast  lag, 
rapidly  at  first,  and  then  approach  constant  values  (approximately  15  cm/s  for  v  and  20  cm/s  for 
u )  at  large  lags.  Forecasts  of  the  northward  component  are  more  accurate  than  those  of  the 
eastward  component  for  all  forecast  lags,  with  the  difference  being  of  the  order  of  3-6  cm/s.  This 
is  likely  a  regional  effect,  resulting  from  the  fact  that  the  wind-driven  currents  are  predominantly 
in  the  east-west  direction,  whereas  the  more  easily  predicted  tidal  currents  are  oriented  more 
north-south.  The  performance  of  the  Gauss-Markov  method  is  clearly  superior  to  the  hedging 
method  at  all  forecast  lags  by  approximately  2-3  cm/s. 


Figure  C-6.  U  RMS  differences  averaged  over  the  entire  domain  as  a  function  of  forecast  lag. 

The  dashed  and  solid  curves  are  the  results  respectively  of  using  the  hedging 
method  and  the  Gauss-Markov  technique  for  estimating  the  non-tidal  current. 
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C.5.2  Spatial  Structure  of  the  Differences 

In  the  FRONT  region  the  dynamics  and  geometry  conspire  to  create  substantial  spatial  variations 
in  the  magnitude  of  the  tidal  and  non-tidal  current  components.  It  is  important,  therefore,  to 
investigate  how  the  algorithm  performs  in  different  parts  of  the  region.  This  is  accomplished  by 
evaluating  the  URMS  at  each  CODAR  grid  point  by  setting  N=1  in  and  t„  =  6,  12  and  24  hours  in 

equation  (4).  Maps  of  the  structure  of  U RMS  are  presented  in  Figure  C-7,  Figure  C-8  and  Figure 
C-9. 
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Figure  C-7.  Spatial  structure  of  the  RMS  differences  for  predictions  made  using  the  hedging 
method  (left)  and  the  Gauss-Markov  method  (right),  evaluated  6  hours  after  the 
prediction  was  made. 
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Figure  C-8.  Spatial  structure  of  the  RMS  differences  for  predictions  made  using  the  hedging 
method  (left)  and  the  Gauss-Markov  method  (right),  evaluated  12  hours  after  the 
prediction  was  made. 
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U  RMS  Diff.  (Hedging),  At=24  hrs 
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Figure  C-9.  Spatial  structure  of  the  RMS  differences  for  predictions  made  using  the  hedging 
method. 


It  is  immediately  apparent  from  Figure  C-7,  Figure  C-8  and  Figure  C-9  that  the  spatial  structures 
in  the  error  maps  for  the  two  methods  are  quite  similar.  As  noted  above,  differences  in 
magnitude  are  observed  at  all  lags  with  the  Gauss-Markov  method  performing  somewhat  better 
for  both  components.  The  largest  errors  in  the  eastward  component  occur  along  the  edges  of  the 
domain,  especially  the  areas  south  of  Montauk  Point  and  Block  Island.  These  are  regions  where 
the  CODAR  measurements  themselves  are  least  accurate  with  respect  to  the  eastward  component 
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(see  Figure  C-4,  top).  As  will  be  discussed  below,  the  increased  measurement  uncertainty 
damages  both  the  prediction  itself  and  its  evaluation.  The  spatial  structure  of  the  error  in  the 
predicted  northward  component  is  quite  different.  Although  predictions  are  poor  at  a  few 
CODAR  grid  points  along  the  edge  of  the  coverage  region  southeast  of  Block  Island,  generally 
the  largest  errors  occur  in  the  zone  between  Montauk  Point  and  Block  Island.  Because  of  the 
coastal  geometry,  the  strongest  north-south  non-tidal  currents  and  vertical  stratification  occur  in 
this  area,  suggesting  a  linkage  between  current  magnitude  and  the  prediction  error. 

C.5.3  Consistency  of  the  Uncertainty  estimates 

For  the  current  predictions  to  be  useful  for  the  prediction  of  Lagrangian  trajectories  and  the 
uncertainty  of  the  target  location,  it  is  essential  that  the  uncertainties  in  the  predicted  current  at 
the  time  that  the  predictions  are  made  be  established.  Because  of  the  superior  performance  of 
the  Gauss-Markov  technique  and  the  fact  that  it  provides  error  estimates  as  well,  we  focus  in  the 
following  on  the  predictions  made  using  this  method. 

The  RMS  differences  between  the  predicted  and  observed  surface  currents  result  from  three 
contributions.  The  observations  themselves  are  uncertain,  and  there  are  errors  associated  with 
the  prediction  of  the  tidal  and  the  non-tidal  motions.  Mathematically,  this  can  be  expressed: 
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RMS 
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®  tidal 


+  <7 


gm  5 


(5) 


where  URMS  is  the  difference  between  the  observations  and  predictions,  and  <rofc,  oiidai,  and  crgm 

are  the  uncertainties  in  the  observations,  the  tidal  predictions,  and  the  Gauss-Markov  predictions 
of  the  non-tidal  current,  respectively.  To  evaluate  the  consistency  of  the  observations  and  the 
assumptions,  we  compare  the  left-hand  side  of  (5)  with  the  right-hand  side  for  both  east  and 
north  components  in  Figure  C-10.  The  fact  that  the  points  generally  fall  along  line  of  unit  slope, 
(the  red  lines  in  Figure  C-10),  is  evidence  that  the  prediction  errors  are  consistent  with  the 
assumptions  and  estimates  of  the  measurement  and  algorithm  errors. 

Equation  5  can  be  rearranged  to  provide  an  expression  for  the  total  prediction  error  in  terms  of 
the  actual  RMS  difference  and  the  estimated  observational  error: 
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a  d  provides  an  a  posteriori  estimate  of  the  true  forecast  performance  since  the  prediction  can 
not  be  better  than  the  observations.  Figure  1 1  shows  the  evolution  of  <rpred  for  the  Gauss- 
Markov  method  averaged  over  the  whole  domain.  The  prediction  error  is  only  about  2  cm/s  less 
than  the  RMS  observation/prediction  difference,  reflecting  the  fact  that  <y2obs  in  Figure  C-4  is 

much  smaller  than  U2RMS  in  Figure  C-7,  Figure  C-8  and  Figure  C-9.  There  is,  therefore, 

considerable  scope  for  improvement  in  the  forecast  model  by  the  inclusion  of  dynamical 
constraints.  The  prediction  errors  for  u  and  v  are  10  cm/s  or  less  at  lags  less  than  2  hours  and 
increase  to  18  cm/s  and  13  cm/s  for  u  and  v  respectively  at  forecast  lag  of  25  hours. 
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Figure  C-10.  Observed  RMS  difference  between  observed  and  predicted  currents  at  12  hours 
versus  the  sum  of  the  errors  in  the  observations  and  in  the  predictions  (tidal  and 
non-tidal)  (a),  for  the  eastward  component  and  (b),  for  the  northward  component. 


Prediction  Error  vs.  Time  Since  Prediction 


Time  Since  Prediction  (h) 


Figure  C-l  1.  The  solid  curves  are  the  root  mean  square  (RMS)  difference  between  predicted  and 
observed  current,  averaged  over  the  entire  CODAR  domain,  as  a  function  of  the 
time  since  the  prediction  using  the  Gauss-Markov  method.  The  dashed  curves  are 
the  estimated  prediction  errors. 
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C.5.4  Sources  of  Error 


The  observations  and  predictions  at  the  two  CODAR  grid  points  shown  in  Figure  1  will  now  be 
studied  to  examine  the  sources  of  error  in  the  predictions.  Figure  C-12  (a)  shows  the  evolution 
of  the  east  component  of  the  current  observed  at  station  1  in  blue  and  the  predicted  current  with  a 
forecast  lag  of  12  hours  (i.e.  the  forecast  computed  at  time  t-2  hours)  in  red.  The  corresponding 
north  components  are  shown  in  (c).  Station  1  was  chosen  to  represent  the  areas  where  tides 
dominate.  In  contrast,  station  2  is  representative  of  areas  in  which  tides  are  weak  and  the  motion 
is  dominated  by  fluctuations  at  longer  time  scales.  Figure  C-12  (e)  and  (g)  compare  the  east  and 
north  components  respectively. 

Time  series  of  differences  between  the  observations  and  predictions  at  station  1  (Figure  C-12b 
and  d)  do  not  exhibit  the  strong  tidal  fluctuations  present  in  the  observations  (Figure  C-12a  and 
c).  This  suggests  that  the  tidal  prediction  is  working  well.  This  is  consistent  with  the  estimated 
tidal  prediction  errors  in  Figure  C-4.  Major  errors  arise  from  the  non-tidal  prediction.  Visual 
comparison  of  the  difference  time  series  at  the  two  stations  (Figure  C-12,  b  and  f,  and  d  and  h), 
indicates  a  correlation,  implying  similar  errors  tend  to  occur  simultaneously  at  the  two  stations 
and  that  there  is  a  common  source. 

Prior  analysis  of  the  CODAR  observations  has  shown  that  the  wind  is  effective  in  driving  surface 
current  fluctuations  over  the  whole  domain  (Ullman  and  Codiga,  2004).  To  examine  this 
mechanism  we  obtained  wind  measurements  from  the  NDBC  Buoy  44017,  located 
approximately  23  miles  southwest  of  Montauk  Point.  Figure  C-13  shows  the  dependence  of  the 
prediction-observation  differences  (with  a  forecast  lag  of  12  hours)  on  the  wind  at  the  time  of  the 
current  observation.  Positive  correlations  between  the  differences  in  the  eastward  currents  and 
the  eastward  wind  velocity,  and  between  the  northward  wind  and  northward  velocity  components 
are  apparent.  For  example,  strong  eastward  winds  tend  to  produce  an  eastward  flow  that  is 
stronger  than  the  predicted  eastward  flow.  This  results  in  a  positive  difference  ( uobs  -  upred ). 

The  behavior  is  similar  for  northward  winds. 

We  performed  a  similar  analysis  using  the  observation-prediction  differences  at  a  forecast  lag  of 
3  hours.  These  display  weaker  correlation  with  the  wind.  This  can  be  understood  by  noting  that 
the  Gauss-Markov  prediction  tends  toward  the  time  series  mean  as  the  forecast  lag  increases.  So 
for  12-hour  forecast  lag,  the  prediction  is  closer  to  the  mean  than  for  the  3 -hour  case,  so  the 
wind-driven  current  fluctuations  will  produce  larger  differences  that  are  correlated  with  the  wind. 
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Figure  C-12.  Comparison  of  CODAR  currents  in  December  2002  with  predicted  currents  using 
the  Gauss-Markov  method  12  hours  prior  to  the  observation  times.  Shown  are 
results  from  stations  1  (a-d)  and  2  (e-h)  in  Figure  C-l.  Panels  (a)  and  (e)  show  the 
observed  and  predicted  u  (east);  (b)  and  (f)  show  the  difference  between  the 
observed  and  predicted  u;  (c)  and  (g)  show  the  observed  and  predicted  v  (north); 
(d)  and  (h)  show  the  difference  between  the  observed  and  predicted  v. 
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Figure  C-13.  Differences  between  observed  and  predicted  (using  the  Gauss-Markov  method  and 
evaluating  at  forecast  lag  of  12  hours)  currents  versus  the  wind  measured  at 
National  Data  Buoy  Center  buoy  44017.  (a,  c)  eastward  currents  at  stations  1  and  2 
(Figure  B-l)  versus  the  eastward  wind,  (b,  d)  northward  currents  at  stations  1  and 
2  versus  the  northward  wind.  Correlation  coefficients  are  given  above  each  plot. 

C.6  Trajectory  Predictions  Using  CODAR  data  and  Current  Predictions 

To  demonstrate  the  potential  practical  utility  of  our  empirical  forecast  approach,  we  now  show 
its  performance  in  a  problem  of  significant  importance.  The  ability  to  predict  the  trajectory  of  an 
object  released  in  the  ocean  at  a  known  time  and  location  is  critical  to  effective  marine  search 
and  rescue  operations.  It  also  has  application  in  the  management  of  oil  spill  and  some  military 
operations.  Since  we  have  developed  a  forecast  capability  for  the  Eulerian  currents,  we  now 
evaluate  its  utility  in  the  prediction  of  particle  motion.  We  first  introduce  the  method  we  employ 
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for  the  Euler-Lagrange  transformation  and  then  test  the  approach  by  comparing  drifter  tracks  to 
trajectory  predictions  based  on  forecast  currents,  and  ones  based  on  current  observations.  These 
comparisons  allow  us  to  estimate  the  expected  errors  and  the  relative  contributions  of  the 
uncertainty  in  the  data  and  the  forecast  algorithm.  We  also  evaluate  the  effectiveness  of  an 
algorithm  like  that  used  at  the  moment  to  quantitatively  establish  the  magnitude  of  the 
improvement  that  could  be  achieved  with  an  SCMS. 

C.6.1  Trajectory  Prediction 

By  definition,  the  Eulerian  and  Lagrangian  velocity  at  any  particular  point,  xp,  and  time,  t,  the 
Lagrangian  velocity  and  Eulerian  velocity,  u(xp,t),  must  be  equal,  i.e., 


dt 


(7) 


The  SCMS  system  provides  us  with  uObs(x,0  at  hourly  intervals  and  the  forecast  algorithm  allows 
us  to  compute  upred(x,t)  for  a  short  time  into  the  future.  We  can  use  either  of  these  to 
approximate  the  right-hand  side  of  (7).  Obviously  the  forecast  must  be  used  operationally,  but 
since  the  predicted  currents  can  not  be  shown  to  be  any  more  accurate  than  the  observations,  use 
of  the  data, 


u=u„bs(x,0,  (8) 

in  (7)  establishes  the  upper  bound  for  the  performance  of  the  prediction  of  the  Lagrangian 
trajectory. 

The  numerical  integration  of  (7)  is  quite  straightforward  and  we  employed  a  simple  second-order 
predictor-corrector  scheme  (Press  et  al,  1992).  The  velocity  must  be  interpolated  from  the 
data/prediction  grid  to  xp  and  we  employed  a  nearest  neighbor  averaging  scheme.  The  velocity 
components  at  the  closest  four  grid  points  are  weighted  inversely  as  the  square  of  the  distance 
from  xp.  This  is  a  robust  method  in  the  sense  that  it  will  provide  an  interpolated  velocity  in  the 
presence  of  data  gaps.  In  the  absence  SCMS  data,  archived  tidal  constituents  at  the  nearest 
station  available,  xN,  could  be  used  to  estimate  the  right  side  of  (7).  This  is  equivalent  to  the 
approximation 


u(x^)=u(w)=uv(0  (9) 

which  ignores  the  spatial  structure  of  the  circulation.  Since  this  is  all  that  is  available  to  the  CG 
at  this  time  we  compare  the  consequences  approximation  (9)  to  our  approach,  u(xp,t)  = 

Upred(xp,0>  and  to  the  best  prediction  possible  with  SCMS  data,  (8).  Note  that  the  tidal 
constituents  used  in  this  study  were  obtained  from  the  National  Oceanic  and  Atmospheric 
Administration  for  the  location  located  41°2.9'  N,  071°  57.6'  W. 

An  example  of  24-hour  trajectories  predicted  using  uN  and  u0bs  in  (7),  together  with  the  track  of 

a  drifter  released  near  xN  at  00:00  on  16  December,  2002  are  shown  in  Figure  C-14.  The  release 
point  is  shown  by  the  black  X.  The  surface  velocity  field  at  the  end  of  the  24  hours  is  also  shown 
to  illustrate  the  spatial  structure  of  the  flow  in  the  region.  The  drifter  trajectory  is  represented  by 
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Figure  C-14.  Sample  Trajectory  Predictions. 


the  black  line,  and  the  simulated  tracks  using  and  u0bs  are  colored  red  and  blue  respectively. 
The  open  circles  denote  the  final  positions.  In  this  example  the  trajectory  predicted  using  u N(t) 

has  a  large  error  after  24  hours.  The  separation  between  the  prediction  and  the  drifter  location  is 
approximately  25  km.  In  contrast,  the  trajectory  predicted  using  the  u0bs(x,/j  is  only  3  km.  This 
case  clearly  illustrates  the  consequences  of  neglecting  the  spatial  structure  of  the  circulation.  Of 
course,  in  the  performance  of  the  trajectory,  simulations  are  sensitive  to  the  release  location  and 
the  Figure  14  is  an  extreme  case.  We  explore  the  statistics  of  the  errors  in  the  next  section. 

C.6.2  Trajectory  Error  Statistics 

As  Section  C.3,  each  drifter  trajectory  was  segmented  into  24-hour  segments  with  the  beginning 
of  each  segment  offset  by  6  hours.  Current  predictions  were  generated  for  24  hours  for  each 
drifter  segment.  Each  segment  was  then  treated  like  an  independent  drifter  release  and  trajectory 
predictions  computed  using  three  estimates  of  the  Eulerian  current,  u N(t) ,  u0bS(x,/),  and 

Upred(xT).  This  approach  allowed  140  release  points  to  be  sampled  throughout  the  FRONT 
CODAR  domain.  At  each  time  step,  the  separation  between  the  actual  and  predicted  SLDMB 
positions  was  computed  and  stored.  The  resulting  database  of  separation  versus  time  was  used  to 
produce  the  histograms  of  separation  at  3,  6,  12  and  24  hours  presented  in  Figure  C-15.  We  then 
computed  the  separation  value  that  was  greater  than  95  percent  of  the  estimates  and  present  it  as 
a  time  series  in  Figure  C-16. 
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Figure  C-15.  Block  Island  Region  Trajectory  Error  Histograms. 
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Bl  Region:  95th  Percentile  Separation 


Figure  C-16.  Block  Island  Region  95th  Percentile  Separation. 


The  left  column  in  Figure  C-15  shows  histograms  of  the  separation  between  drifters  and  the 
simulations  (the  target  error)  using  observed  currents,  uobs(x,0-  The  center  and  right  columns 
show  the  performance  of  the  simulations  based  in  the  forecast  currents  upred(x,/),  and  un(0 
respectively.  The  rows  show  results  at  prediction  times  of  3,  6,  12,  and  24  hours. 

Comparison  of  the  histograms  at  any  time  shows  that  the  lowest  modal  values  of  the  target  errors 
occur  when  observed  velocities  are  used  (left  column).  The  central  and  right  columns  show  that 
the  modal  target  error  and  width  of  the  histogram  are  both  substantially  smaller  when  the 
forecast  algorithm  is  employed.  For  example,  at  an  elapsed  time  of  3  hours,  the  upred(x,f) 
predicted  positions  are  always  within  10  km  of  the  actual  drifter  position,  whereas  a  significant 
number  of  the  un(0  based  predictions  are  in  error  by  10-20  km. 

The  evolution  of  the  target  error  statistics  is  most  clearly  summarized  in  Figure  16  by  the  time 
series  of  the  95  percentiles.  The  three  curves  represent  the  three  current  sources.  The  green, 
blue  and  red  lines  show  the  values  obtained  with  u0bs(x,t),  upred(x,0,  and  un(0-  At  all  times,  the 
predictions  made  using  data  are  the  best  and  those  made  using  un(0  are  the  worst,  with 
trajectories  predicted  using  the  forecast  currents  intermediate  in  accuracy.  The  95th  percentile 
separation  value  increases  monotonically  for  the  u0bs(x,/),  and  upred(xT)  based  trajectory 
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predictions  with  values  at  24  hours  of  15  and  30  km  respectively.  Using  Un(/),  the  separation 
decreases  slightly  between  about  6  and  12  hours  before  rising  again  to  a  24-hour  value  of  40-45 
km.  The  reason  for  this  decrease  at  intermediate  times  is  due  to  the  fact  that  the  uN(0  currents 
are  mainly  tidal  with  a  period  of  12.4  hours.  If  applied  in  a  region  where  tidal  currents  are 
actually  weak  (station  2  in  Figure  C-l  for  example),  the  predicted  drifter  motion  will  show  a 
large  oscillation  with  this  period,  whereby  the  actual  drifter  will  move  slowly  more  or  less  in  one 
direction.  At  the  end  of  12.4  hours  in  this  case,  the  simulated  drifter  will  have  almost  returned  to 
its  starting  point  and  therefore  be  closer  to  the  actual  drifter  position. 

C.6.3  Trajectory  Errors  versus  Wind 

CODE-type  drifters  are  thought  to  exhibit  negligible  windage  (motion  relative  to  the  surface 
water  driven  because  of  wind).  This  was  investigated  empirically  by  comparing  drifter 
prediction  errors,  in  the  Block  Island  region,  with  winds  measured  at  a  nearby  NOAA  buoy 
(44017).  The  separation  between  predicted  and  actual  drifter  position  after  24  hours  was  plotted 
versus  the  mean  wind  magnitude  over  the  24-hour  prediction  interval,  as  shown  in  Figure  C-l 7. 
The  three  panels  show,  from  top  to  bottom,  the  comparisons  of  wind  speed  to  separation  for 
predictions  using  u0bs(x,t),  upred (x,t),  and  un(0-  The  first  comparison  shows  that  prediction 
errors  using  observed  currents,  u0bs(x,f),  for  the  prediction  do  not  correlate  with  the  wind.  We 
interpret  that  this  is  because  the  CODAR  observations  include  the  real  wind-driven  surface 
current  component,  which  is  known  from  prior  unpublished  work  to  be  large.  In  contrast,  the 
drifter  trajectories  simulated  by  the  predicted  currents,  upred(x,t)  and  un (t),  tend  to  be  more 
seriously  in  error  under  high  wind  conditions.  Since  uN(t),  does  not  include  any  subtidal 
frequency  motion,  it  is  not  surprising  that  the  target  errors  correlate  with  the  wind  speed. 
However,  the  CODAR  observation,  and  the  predictions  based  on  them,  upred(x,t),  do  have 
information  about  the  wind  driven  current,  though  it  appears  that  the  algorithm  could  be 
improved  by  capturing  the  effect  of  wind  more  directly.  Note  however,  that  the  drifters  were  in 
the  region  of  FRONT  CODAR  coverage  for  only  10  days  and  we  may  not  have  adequately 
sampled  the  response  to  wind. 

C.7  Conclusions  and  Recommendations 

The  surface  current  forecast  algorithm  developed  in  this  paper  is  based  on  the  independent 
prediction  of  the  tidal  and  non-tidal  surface  current  components.  Tidal  currents  are  predicted 
using  harmonic  constituents  derived  from  a  1 -month  record  of  CODAR  estimated  surface 
currents.  Two  methods  for  predicting  the  non-tidal  component  were  tested:  hedging  and  Gauss- 
Markov  estimation.  The  Gauss-Markov  method  was  clearly  superior  at  all  forecast  lags  up  to  25 
hours,  the  maximum  lag  that  was  tested.  Since  the  calculation  of  the  autocorrelation  function  for 
the  Gauss-Markov  method  requires  at  least  a  30-day  archive  of  currents,  the  hedging  approach 
may  be  useful  at  new  SCMS  locations.  Either  of  these  strategies  is  capable  of  turning  data  from 
a  coastal  current  observing  network  into  useful  forecast  products. 

Analysis  of  the  combined,  tidal  plus  Gauss-Markov  non-tidal,  prediction  algorithm  over  a  1- 
month  period  gives  estimates  of  the  region-averaged  prediction  error  of  9-18  cm/s  for  the 
eastward  component  and  7-13  cm/s  for  the  northward  component.  The  range  of  prediction  errors 
reflects  the  increase  in  prediction  error  with  forecast  lag.  Note  that  the  performance  should  be 
expected  to  be  site-specific. 
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Figure  C-17.  Block  Island  Region  24-hour  Separation  versus  Wind.  A  key  component  of  the 
Short  Term  Predictive  System  (STPS)  is  the  estimation,  at  the  time  of 
prediction,  of  the  uncertainties  associated  with  the  forecast  velocities.  The 
magnitudes  of  the  predicted  uncertainties  are  consistent  with  the  computed 
RMS  differences  and  the  estimated  uncertainties  in  the  observations  themselves. 
This  leads  to  the  conclusion  that  the  forecast  uncertainties  can  be  reliably  used 
to  characterize  the  statistics  of  errors  in  the  predictions.  The  prediction 
uncertainties  are  essential  to  modeling  Lagrangian  dispersion. 

The  correlation  of  observation/prediction  differences  with  wind,  which  is  known  to  be  the  major 
forcing  function  for  non-tidal  fluctuations  in  the  FRONT  region,  is  somewhat  surprising.  It 
suggests  that  the  relative  weight  placed  on  the  data  contributing  to  the  prediction  is  not  quite 
correct.  This  could  result  from  incorrect  characterization  of  either  the  observational  errors  in  the 
velocity  components  or  the  data  autocovariance  functions.  The  observational  errors  are 
estimated  by  assuming  that  the  uncertainty  in  the  radial  velocity  from  a  given  site  is  spatially  and 
temporally  constant,  which  has  yet  to  be  verified  for  the  CODAR  system.  The  autocovariance 
functions  used  in  our  implementation  of  the  Gauss-Markov  estimation  procedure  are  averaged 
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spatially  over  the  CODAR  domain  and  are  then  scaled  by  the  actual  data  variance  at  each  grid 
point.  This  spatial  averaging  is  done  to  produce  smoother  autocovariance  functions,  the  use  of 
which  helps  to  reduce  the  incidence  of  wild  predictions  at  grid  points  with  frequent  data  gaps. 
However,  there  is  clearly  some  spatial  structure  present  in  the  autocovariance  functions  and  its 
neglect  can  lead  to  misrepresentation  of  the  data  covariance  in  some  locations  in  the  region. 
Further  work  is  necessary  to  determine  whether  the  use  of  more  robust  autocovariance  functions 
derived  from  historical  data  records  of  say  several  months  will  produce  more  accurate 
predictions. 

Although  more  work  is  necessary  to  confirm  this  point,  preliminary  work  on  the  prediction 
methodology  carried  out  using  CODAR  data  from  August  2001  resulted  in  somewhat  lower 
prediction  error  during  that  time  period  as  compared  with  the  December  2002  period  presented  in 
this  report.  As  the  accuracy  of  the  tidal  prediction  is  unlikely  to  vary  in  a  systematic  manner 
over  the  seasons,  this  suggests  that  the  non-tidal  predictability  may  be  reduced  in  winter.  This 
may  be  due  to  increased  non-tidal  (wind-driven)  current  fluctuations  driven  by  wintertime  wind 
intensification. 

Comparison  of  drifter  trajectory  predictions  with  observed  trajectories  demonstrates  that  the 
SCMS-based  forecasts  are  clearly  superior  to  current  practice.  This  arises  from  two  major 
improvements:  a  good  representation  of  the  spatial  structure  of  the  mean  and  tidal  flow,  and 
inclusion  of  an  estimate  of  the  low  frequency  variability.  The  simple  Gauss-Markov  approach 
results  in  mean  drifter  position  errors  after  6  hours  of  7  km  and  13  km  after  12  hours.  These 
errors  are  approximately  half  of  that  obtained  using  the  tidal  current  predictions  from  the 
archives  maintained  by  NOAA,  as  is  currently  the  operational  practice  in  the  CG.  Operational 
use  of  the  currently  available  technology  could,  therefore,  have  the  potential  of  reducing 
predicted  search  areas  by  a  factor  of  four. 
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Abstract 

An  important  aspect  of  trajectory  modeling  in  support  of  search  and  rescue  applications  is  the 
assessment  of  the  uncertainty  in  the  search  target  position.  Many  factors  contribute  to  the 
uncertainty  and  spread  of  the  final  position  distribution,  including  uncertainty  in:  initial  position, 
initial  time  of  incident,  wind  and  leeway,  and  surface  currents.  This  paper  focuses  on  the 
contribution  of  the  uncertainty  in  the  surface  currents  to  the  search  and  rescue  final  probability 
distribution.  Monte-Carlo  particle  trajectory  simulations  using  surface  currents  derived  from 
CODAR  High  Frequency  (HF)  radar  systems  were  performed  with  sub-gridscale  motions 
parameterized  using  random  walk  and  random  flight  turbulence  models,  and  the  resulting 
trajectories  were  compared  to  the  trajectories  of  surface  drifting  buoys.  Velocity  statistics  for  the 
turbulence  models  were  derived  from  covariance  functions  of  time  series  of  differences  between 
CODAR  and  drifter  estimates  of  surface  currents.  Drifter  positions,  predicted  over  24-hour 
periods  using  CODAR  currents,  were  found  to  be  more  accurate  predictors  of  the  real  drifter 
positions  than  were  initial  positions.  The  particle  dispersion  predicted  by  24-hour  Monte-Carlo 
simulations  using  the  estimated  fluctuation  statistics  was  used  to  determine  95  percent 
confidence  regions,  which  were  subsequently  tested  using  the  real  drifter  location.  Search  areas 
determined  from  the  95  percent  regions  using  the  random  flight  model  were  found  to  be  80-90 
percent  effective,  in  reasonable  agreement  with  the  expected  95  percent  success  rate.  On  the 
other  hand,  search  areas  computed  from  random  walk  modeling  were  found  to  be  inadequate 
unless  the  diffusion  coefficient  was  increased  to  approximately  the  random  flight  value.  The 
effectiveness  of  the  random  flight  search  areas  lends  support  to  the  proposed  methodology  for 
estimating  turbulence  parameters  from  drifter-CODAR  velocity  differences. 

D.l  Introduction 

The  recent  proliferation  of  coastal  High  Frequency  (HF)  surface  wave  radar  installations  for 
mapping  surface  ocean  currents  provides  a  rapidly  expanding  capability  for  real-time  observation 
of  surface  currents.  These  data  have  the  potential  to  dramatically  improve  the  efficiency  and 
success  rate  of  search  and  rescue  operations  in  coastal  waters.  For  this  reason,  the  R&DC  has 
initiated  a  program  to  assess  the  effectiveness  of  trajectory  predictions  using  currents  derived 


D-l 


from  HF  radar,  and  if  warranted,  to  implement  the  use  of  this  technology  on  an  operational  basis. 
Preliminary  work  toward  this  end  was  reported  by  Ullman,  O’Donnell,  Edwards,  Fake, 
Morschauser,  Sprague,  Allen,  Krenzien,  (2003). 

Trajectory  modeling  applied  for  search  and  rescue  operations  must  provide  a  measure  of  the 
uncertainty  in  the  surface  current  portion  of  drift  of  the  search  and  rescue  object  in  order  for  an 
optimal  search  area  to  be  delineated.  The  identification  of  a  search  area  is  based  on  a  balance 
between  the  need  to  define  a  large  enough  area  to  ensure  that  the  target  is  within  the  search  area 
and  the  reality  that  search  resources  are  finite.  The  size  of  an  operational  search  area  is  related  to 
the  magnitude  of  various  uncertainties  including  those  introduced  by  poorly  known  initial  target 
position  and  time,  velocity  errors,  and  target  leeway.  This  paper  focuses  on  the  estimation  of  a 
search  area  due  only  to  uncertainties  in  the  drifter  advective  velocity  provided  by  HF  radar 
systems. 

Search  areas  are  estimated  in  this  paper  using  a  Monte-Carlo  approach  whereby  the  trajectories 
of  an  ensemble  of  1,000  particles  are  modeled  over  24  hours.  The  advective  velocity  is 
decomposed  into  a  deterministic,  large-scale  component  measured  by  the  CODAR  plus  a  non- 
deterministic  component  representing  a  combination  of  sub-gridscale  motion  and  errors  in  the 
CODAR  velocity.  We  make  use  of  the  hierarchy  of  stochastic  particle  models  outlined  by  Griffa 
(1996),  testing  the  first  two  of  these,  the  random  walk  and  random  flight  models,  in  the  present 
work. 

D.2  Data  and  Methodology 

D.2.1  Drifter  Trajectories 

A  number  of  drifter  releases  in  the  Mid  Atlantic  Bight  were  performed  by  the  R&DC  over  the 
period  2002-2004  with  the  objective  of  providing  a  data  set  with  which  to  assess  trajectory 
predictions.  Drifters  were  released  in  Block  Island  Sound  within  the  coverage  region  of  the 
standard-range  CODAR  system  operated  by  the  Universities  of  Rhode  Island  and  Connecticut 
(Figure  D-l).  A  number  of  these  drifters  eventually  passed  through  the  coverage  region  of 
Rutgers  University’s  long-range  CODAR  system  and  a  number  of  additional  drifters  were 
released  within  that  zone  as  well  (Figure  D-2).  The  CODE-type  drifters  (Allen,  1996)  recorded 
Global  Positioning  System  (GPS)  fixes  on  one-half  hour  intervals  and  transmitted  positions  via 
Argos.  With  one-half  hour  sampling  and  a  nominal  GPS  position  uncertainty  of  10  m,  the 
uncertainty  in  the  velocity  of  the  drifter  is  0(1  cm/s),  which  will  be  seen  to  be  small  compared  to 
the  observed  differences  between  drifter  and  radar-derived  velocities. 
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Figure  D-l.  Trajectories  of  surface  drifters  deployed  during  December  2002  and  March  2003  in 
the  Block  Island  CODAR  region.  The  black  dots  show  the  release  points  of  each 
drifter  (note  that  some  were  retrieved  and  redeployed).  The  green  diamonds  show 
the  locations  of  the  CODAR  sites,  and  the  dashed  line  shows  the  approximate  10 
percent  coverage  region. 
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Figure  D-2.  Trajectories  of  surface  drifters  deployed  during  March  2003  (blue)  and  July  2004 
(red)  that  passed  through  the  New  Jersey  shelf  CODAR  domain.  The  black  dots 
show  the  release  points  of  each  drifter.  The  green  diamonds  show  the  locations  of 
the  CODAR  sites.  The  dotted  line  shows  the  nominal  coverage  zone  and  the  dashed 
line  shows  the  approximate  10  percent  coverage  region. 


Drifters  were  deployed  in  December  2002  in  the  Block  Island  region  and  again  in  March  2003. 
A  number  of  the  latter  group  subsequently  moved  southwest  on  the  shelf  and  passed  through  the 
New  Jersey  shelf  CODAR  coverage  region.  A  deployment  of  drifters  was  also  made  in  March 

2003  within  the  New  Jersey  shelf  domain.  A  final  set  of  deployments  was  made  during  July 

2004  in  both  regions;  however  the  drifters  released  in  the  Block  Island  region  rapidly  exited  that 
domain  and  did  not  provide  significant  trajectory  segments.  Table  D-l  summarizes  the  drifter 
deployments.  Initial  analysis  focused  on  drifters  within  the  Block  Island  region  during  2002- 
2003  and  on  those  within  the  New  Jersey  shelf  region  during  2004.  The  New  Jersey  shelf 
drifters  from  2003  were  subsequently  used  for  the  purpose  of  verifying  the  choice  of  turbulence 
parameters  made  using  the  2004  data. 
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Table  D-l .  R&DC  Mid  Atlantic  Bight  drifter  releases  during  2002-2004. 


Date 

Region 

Number 

Notes 

Dec.  16-18,  2002 

Block  Island 

9 

Retrieved  and  redeployed  several. 

Mar.  27,  2003 

Block  Island 

4 

Passed  through  NJ  Shelf  region. 

Mar.  27,  2003 

NJ  Shelf 

8 

Jul.  27,  2004 

Block  Island 

3 

Rapidly  left  domain. 

Passed  through  NJ  Shelf  region. 

Jul.  27,  2004 

NJ  Shelf 

4 

D.2.2  CODAR  Surface  Currents 

HF  radar  surface  currents  were  obtained  with  CODAR  SeaSonde  systems  located  in  the  Mid 
Atlantic  Bight  in  the  region  around  Block  Island  and  on  the  shelf  east  of  New  Jersey.  A  three- 
site  standard-range  (25  MHz)  system  in  the  Block  Island  region  provided  hourly  surface  currents 
at  1.5  km  resolution  over  the  region  shown  in  Figure  D-l.  Radial  velocities  observed  at  each  of 
the  three  sites  were  obtained  with  the  Multiple  Signal  Classification  (MUSIC)  algorithm 
(Schmidt,  1986)  using  measured  antenna  patterns;  and  combining  into  vector  velocities  was 
performed  using  the  manufacturer’s  (CODAR  Ocean  Sensors)  software.  The  New  Jersey  shelf 
surface  currents  with  spatial  resolution  of  6  km  were  provided  by  a  four-site  long-range  (5  MHz) 
system  in  the  area  shown  in  Figure  D-2.  Radial  velocities,  averaged  over  3 -hour  periods  using 
measured  antenna  patterns,  were  output  on  an  hourly  basis  and  were  combined  using  the  Naval 
Postgraduate  School’s  HF  Radar  Toolbox.  A  screening  methodology,  utilizing  a  threshold  on 
the  estimated  geometric  dilution  of  precision  (GDOP)  (Gurgel,  1994),  was  used  to  remove 
current  vectors  derived  from  combinations  of  radials  from  sites  with  poor  geometry.  Vectors  for 
which  the  GDOP  was  greater  than  an  empirical  threshold  of  1.25  were  eliminated  from 
consideration. 

Estimates  of  the  uncertainties  associated  with  surface  currents  derived  from  CODAR  systems 
have  been  provided  using  in  situ  velocity  observations  from  Acoustic  Doppler  Current  Profilers 
(ADCPs).  Chapman  and  Graber  (1997)  cite  differences  of  0(15  cm/s)  between  HF  radar  current 
estimates  and  in  situ  current  measurements.  However,  as  Kohut  et  al.  (2004)  point  out,  these  error 
estimates  include  a  large  component  that  is  due  to  the  different  spatial  scales  and  depths  sampled 
by  HF  radar  and  ADCPs.  They  estimate  the  intrinsic  CODAR  radial  uncertainty  to  be  of  0(5 
cm/s)  for  the  New  Jersey  shelf  long-range  systems  used  in  the  present  study.  Vector  uncertainties 
for  the  Block  Island  region  were  estimated  to  be  3-15  cm/s  with  the  larger  values  observed  along 
the  outer  boundaries  of  the  coverage  regions  where  the  combining  geometry  is  non-optimal  and 
where  signal  to  noise  ratios  increase  (Ullman  and  Codiga,  2004). 

Interpolation  of  velocities  from  the  CODAR  grid  to  the  location  of  a  drifter  was  performed  using 
a  weighted,  nearest-neighbor  scheme  in  which  velocities  from  the  four  nearest-neighbor  grid 
points  were  weighted  inversely  with  distance.  This  method  is  more  robust  than  bi-linear 
interpolation  because  it  is  more  forgiving  of  missing  CODAR  vectors  and  allows  for 
extrapolation  beyond  the  instantaneous  zone  of  CODAR  coverage.  Although  velocities  could  be 
extrapolated  to  any  distance  beyond  the  coverage  zone,  in  practice,  we  restrict  the  trajectory 
prediction  to  those  trajectory  segments  that  start  within  the  nominal  coverage  zones  shown  in 
Figures  D-l  and  D-2.  Also  shown  in  Figure  D-2  is  the  10  percent  coverage  zone.  This  is  the 
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region  within  which  CODAR  vectors  are  obtained  at  least  10  percent  of  the  time,  computed  over 
the  drifter  deployment  period.  Note  that  in  the  Block  Island  region,  the  10  percent  coverage 
region  and  the  nominal  coverage  region  were  the  same.  The  10  percent  coverage  regions  were 
used  to  perform  more  stringent  screening  of  the  trajectories  to  regions  of  more  reliable  CODAR 
currents. 

D.2.3  Monte-Carlo  Trajectory  Prediction 

The  motion  of  a  particle  in  a  two-dimensional  velocity  field  can  be  described  by  the  equation: 

dt*  .  .  ... 

—  =  u{t,r) ,  (1) 

dt 

where  r=(x,y)  denotes  the  position  of  the  particle  and  u=(u,v)  is  the  Eulerian  velocity  at  position 
r  and  time  t.  The  velocity  can  be  decomposed  into  a  large-scale,  slowly  varying  component,  U 
and  a  component,  ut  representing  sub-gridscale  deviations,  which  will  be  referred  to  as 
turbulence: 


u  =  U  +  ut.  (2) 

Surface  current  mapping  radars  such  as  CODAR  can  provide  estimates  of  V  at  spatial  scales  of 
approximately  1.5  km  (6  km)  and  temporal  scales  of  1  hr  (3  hr)  for  standard-range  (long-range) 
systems,  thus  the  turbulent  component  represents  velocity  fluctuations  on  scales  smaller  than 
these.  The  radar-derived  velocity  is  subject  to  significant  uncertainties  such  that  the  large-scale 
component  can  be  expressed  as: 

V  =  URadr+Su,  (3) 

where  Uradar  is  the  radar  measurement  and  5u  is  the  measurement  error.  Combining  (2)  and  (3), 
the  total  Eulerian  velocity  can  be  written: 

»  =  U radar  +  ",  +  5u  =  U radar  +  >  (4) 

where  «'  includes  both  the  turbulent  velocity  and  the  measurement  error. 

Prediction  of  particle  trajectories  in  a  region  of  HF  radar  coverage  was  achieved  by  integrating 
(1)  using  a  predictor-corrector  scheme  with  the  velocity  given  by  (4).  A  time-step  of  1  hr  was 
used  and  radar  velocities  were  interpolated  to  the  particle  location  using  the  weighted  nearest- 
neighbor  approach.  The  methodology  for  specifying  m'  is  described  in  the  following  section. 

D.2.4  Sub-gridscale  model 

Although  «'  is  a  combination  of  radar  measurement  errors,  unresolved  motion,  and  true 
geophysical  turbulent  fluctuations,  we  hypothesize  that  its  properties  can  still  be  described  by 
models  of  turbulence.  This  is  partly  motivated  by  the  fact  that  there  is  no  easy  way  to  separate 
the  two  components  when  comparing  with  velocity  estimates  from  drifters.  The  fact  that  surface 
drifters  are  routinely  launched  in  support  of  search  and  rescue  operations  motivates  the 
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development  of  a  methodology  to  estimate  the  statistics  of  u'  from  drifter-derived  velocities.  The 
fluctuating  velocity,  «'  will  be  referred  to  below  as  turbulence  with  the  understanding  that  this 
quantity  also  contains  a  significant  measurement  error  component. 

Two  models  of  turbulence  were  examined,  both  of  which  are  members  of  the  hierarchy  of 
stochastic  particle  models  reviewed  by  Griffa  (1996).  The  so-called  random  walk  and  random 
flight  models  both  assume  that  the  particle  position  is  Markovian  and  differ  in  that  the  latter  also 
assumes  that  the  particle  velocity  is  a  Markov  variable.  Physically,  the  difference  is  that  the 
turbulence  has  a  finite  temporal  scale  and  therefore  temporal  correlation  in  the  random  flight 
case  whereas  the  random  walk  model  assumes  that  the  turbulent  time  scale  is  infinitesimal.  For 
equal  turbulent  velocity  variance,  one  therefore  expects  greater  dispersion  in  the  former  case. 


The  two  components  of  the  vector  velocity  u'  are  assumed  independent  in  both  models.  The 
random  walk  formulation  for  u' can  be  expressed  as: 


j’X 

u'  =  <7„ - dw, 

"  dt 


(5) 


where  ou  is  the  velocity  standard  deviation  and  dw  is  a  normally  distributed  random  increment 
with  zero  mean  and  second  moment  (dw  •  dw)  =  2  •  dt  with  dt  the  time  step  for  the  integration  of 

(1).  The  turbulent  time  scale  Tu  in  the  discrete  problem  is  not  actually  infinitesimal  but  is 
constrained  to  equal  dtl 2  to  obtain  consistency  of  (5)  and  the  definition  of  velocity  variance 
(Griffa,  1996). 

The  evolution  of  the  turbulent  velocity  in  the  random  flight  turbulence  model  is  described  by: 
du'  =  -—dt  +  ^-dw.  (6) 

T  ’tYi  y  ' 

u  ±u 


The  first  term  on  the  right  of  (6)  represents  the  “memory”  of  the  turbulence,  decaying  with  time 
scale  Tu  which  is  arbitrary  in  this  case.  The  second  term  is  a  Gaussian  random  impulse  as  in  the 
random  walk  case.  The  autocorrelation  function  of  the  turbulent  velocity  decays  exponentially 
with  an  e-folding  time  or  integral  time  scale  equal  to  Tu  (Griffa,  1996). 

It  is  important  to  note  that  the  diffusion  coefficient  for  particles  in  homogeneous  turbulence  at 
times  large  compared  with  Tu  is  defined  as: 


K 


x 


=  af-T. 


(7) 


(Csanady,  1973).  Since  for  the  random  walk  case,  Tu=dtl 2,  the  dispersion  in  numerical 
implementation  of  this  model,  for  a  given  a2 ,  is  dependent  on  the  time  step  employed. 

D.2.5  Estimating  turbulent  velocity  statistics 

Practical  implementation  of  the  aforementioned  turbulence  models  to  determine  the  random 
velocity  components  in  (4)  requires  specification  of  the  velocity  variance,  and  for  the  random 
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flight  case  the  turbulent  time  scale  Tu.  One  could  prescribe  the  diffusion  coefficient  in  both 
cases,  but  we  argue  that  with  relatively  few  drifters  deployed  and  lacking  large  clusters,  it  is 
more  appropriate  to  estimate  the  fluctuating  velocity  statistics  than  the  diffusion  coefficient. 

From  (4),  the  fluctuating  velocity  is:  u  =u-  Umdar ,  and  we  consider  the  drifter  velocity  to 
approximate  the  true  velocity  u. 

Time  series  of  drifter  positions  at  1/2-hourly  intervals  were  used  with  a  central  difference  scheme 
to  compute  drifter  velocities  at  times  (hourly)  corresponding  to  CODAR  observations.  CODAR 
velocities  were  then  spatially  interpolated  to  the  drifter  location  using  a  weighted  nearest- 
neighbor  method.  Velocities  at  locations  outside  of  the  10  percent  coverage  zones  were 
eliminated  from  consideration.  Time  series  of  u '  for  each  drifter  were  then  computed  by 
differencing  the  drifter  and  CODAR  velocities.  These  time  series  were  used  to  estimate 
covariance  functions,  averaged  over  all  drifters,  for  the  u'  and  v'  components  (Figure  D-3). 


Figure  D-3.  Lagged  covariance  functions  of  time  series  of  differences  between  drifter  velocity 
and  CODAR  velocity  averaged  over  all  drifters  within  the  Block  Island  region 
(left)  and  the  New  Jersey  shelf  region  (right).  The  autocovariances  of  the 
eastward  (u)  and  northward  (v)  components  are  shown  by  the  blue  solid  and  the 
red  dashed  lines,  respectively.  The  green  dash-dot  line  denotes  the  cross¬ 
covariance. 


The  autocovariance  functions  exhibit  rapid  decay  at  lags  of  several  hours  and  display  lower 
amplitude  periodicities  at  the  semi-diurnal  period  in  the  Block  Island  region  and  at 
approximately  the  inertial  period  in  the  New  Jersey  shelf  region  (Figure  D-3).  The  cross¬ 
covariance  is  generally  low,  consistent  with  statistical  independence  between  the  two 
components.  Turbulent  variances  were  estimated  as  the  zero-lag  values  of  the  autocovariance 
functions.  The  turbulent  integral  time  scales  were  estimated  as  one  half  the  time  lag  of  the  first 
zero  crossing  of  the  autocovariance  functions.  The  estimates  of  these  parameters  are  given  in 
Table  D-2. 
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Table  D-2.  Estimates  of  turbulence  parameters  from  autocovariance  functions  of  drifter- 


CODAR  velocity  differences  for  the  two  CODAR  regions. 


Region 

ctu  (cm/s) 

Tu  (h) 

ctv  (cm/s) 

Tv  (h) 

Block  Island 

14.4 

1.5 

11.8 

1.5 

NJ  Shelf 

10.6 

3 

12.3 

3 

D.3  Results 

To  evaluate  the  consistency  of  our  approach  to  modeling  dispersion  and  the  estimates  of  the 
parameters,  we  compared  observed  drifter  trajectories  over  24  hours  to  those  predicted 
retrospectively  using  CODAR  velocities.  The  trajectory  of  each  drifter  was  broken  down  into 
24-hour  segments  with  each  segment  overlapping  the  previous  one  by  12  hours  for  drifters 
within  the  Block  Island  region.  An  example  simulation  is  shown  in  Figure  D-4  where  we  show 
the  95th  percentile  region  for  the  final  psuedo-drifter  location  as  a  gray  polygon.  This  region  was 
estimated  by  computing  the  two-dimensional  histogram  of  the  pseudo-drifter  location,  sorting  the 
spatial  bins  based  on  the  number  within  and  then  summing  over  bins  until  95  percent  of  the  total 
(1000)  was  obtained.  Trajectories  were  only  computed  for  segments  in  which  the  starting 
position  of  the  drifter  was  within  the  CODAR  coverage  region.  In  the  analysis  that  follows,  we 
present  statistics  for  all  computed  trajectories  as  well  as  statistics  of  a  screened  subset  of 


Figure  D-4.  Example  drifter  trajectory  within  the  Block  Island  CODAR  region  showing  the  real 
drifter  path  over  24  hrs  in  red  and  the  CODAR  predicted  position  in  green,  with 
final  positions  denoted  by  the  circles.  The  ensemble  of  trajectory  prediction 
endpoints  using  the  random  walk  (left)  and  random  flight  (right)  models  are  shown 
as  blue  dots.  The  gray  polygons  denote  the  region  within  which  95  percent  of  the 
drifter  final  positions  lie. 
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trajectory  segments  for  which  the  end  position  is  also  located  within  the  coverage  region.  The 
number  of  comparisons  versus  the  time  since  the  start  of  the  prediction  is  shown  in  Figure  D-5. 
The  decrease  in  number  with  time  in  the  Block  Island  region  results  from  the  retrieval  and 
redeployment  of  drifters  that  left  the  CODAR  domain;  thus  a  number  of  “short”  trajectories  (< 
24  hr)  is  present  in  the  trajectory  ensemble. 


Figure  D-5.  The  number  of  comparisons  between  predicted  and  true  drifter  position  versus  time 
for  no  screening  (solid  lines)  and  for  screening  using  the  10  percent  coverage  zone 
(dashed  lines)  for  the  Block  Island  region  (left)  and  New  Jersey  Shelf  region  (right). 


D.3.1  Accuracy  of  Predictions 

The  accuracy  of  a  drifter  trajectory  prediction  is  measured  by  the  distance  between  the  real 
drifter  and  the  pseudo-drifter.  This  was  computed  for  each  hour  of  each  24-hour  trajectory 
segment.  The  ensemble  mean  separation  and  the  95th  percentile  separation  is  presented  in 
Figure  D-6  for  all  trajectories  and  for  the  screened  subset.  Mean  separations  in  the  two  regions 
increase  with  time  in  a  linear  fashion  with  some  indication  that  separations  at  short  times 
increase  at  a  slightly  faster  rate.  After  24  hours,  mean  separation  is  approximately  7  km  (6  km) 
for  the  Block  Island  unscreened  (screened)  ensembles  and  about  8  km  (7  km)  for  the  drifters  in 
the  New  Jersey  shelf  region.  The  95th  percentile  separation  values  are  also  slightly  higher  in  the 
latter  domain,  for  the  unscreened  subset,  reaching  20  km  at  the  end  of  24  hours  compared  to 
about  18  km  for  the  Block  Island  drifters.  When  drifters  leaving  the  CODAR  region  are 
eliminated,  the  95th  percentile  separations  (24  hr)  are  about  12  km  and  15  km  in  the  two  regions 
respectively.  The  general  decrease  in  error  with  screening  is  not  surprising,  and  is  consistent 
with  the  occurrence  of  relatively  large  trajectory  prediction  error  along  the  outer  boundary  of  the 
CODAR  domain  where  radar  velocity  errors  increase. 
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Figure  D-6.  Separation  between  actual  and  predicted  drifter  position  as  a  function  of  time  since 
start  of  prediction,  averaged  over  all  trajectory  segments  that  start  within  the 
nominal  coverage  zone  for  the  Block  Island  region  (upper  left)  and  the  New  Jersey 
shelf  region  (upper  right).  The  blue  curves  show  the  mean  separation  (solid  line 
and  plusses)  and  the  95th  percentile  separation  (solid  line)  of  the  real  drifter  and  the 
predicted  position.  The  red  curves  show  the  mean  distance  (dashed  line  and 
plusses)  and  the  95th  percentile  distance  (dashed  line)  that  the  real  drifter  moved 
over  the  prediction  time.  The  lower  plots  show  the  same  statistical  measures 
averaged  over  all  segments  that  both  start  and  end  within  the  10  percent  coverage 
zone  for  the  Block  Island  region  (lower  left)  and  New  Jersey  shelf  (lower  right). 
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To  place  the  CODAR-based  trajectory  accuracies  into  context,  and  to  evaluate  whether  CODAR 
predictions  provide  an  improvement  over  the  so-called  persistence  forecast  (last  known  position), 
we  also  show,  in  Figure  D-6,  the  mean  and  95th  percentile  values  of  the  distance  traveled  by  the 
real  drifters.  Drifters  released  in  the  Block  Island  region  tend  to  travel  farther  than  those  released 
on  the  open  shelf.  For  the  unscreened  case,  the  mean  (95th  percentile)  distance  traveled  after  24 
hours  in  the  Block  Island  region  is  15  km  (38  km)  compared  with  about  1 1  km  (23  km)  in  the 
New  Jersey  shelf  region.  The  net  result  is  that  drifter  locations  predicted  using  CODAR  currents 
in  the  Block  Island  region  are  on  average  approximately  50  percent  closer  to  the  real  drifter 
position  than  are  the  last  known  positions.  This  effect  is  less  pronounced  on  the  New  Jersey 
shelf  where  the  mean  distance  from  the  predicted  position  to  the  real  drifter  is  about  70  percent 
of  the  distance  traveled.  Similar  conclusions  are  reached  using  the  screened  drifter  subsets, 
although  the  screening  can  be  seen  to  sharply  reduce  the  95th  percentile  separation  value  in  the 
Block  Island  region. 

D.3.2  Uncertainty  Bounds  for  Predictions 

Monte-Carlo  simulation  of  drifter  trajectories  provides  an  ensemble  of  final  drifter  locations  that 
can  be  used  to  construct  a  search  area.  The  search  area  was  defined  as  the  region  in  which  95 
percent  of  the  pseudo-drifters  were  located,  estimated  as  discussed  above.  At  each  hour  of  each 
modeled  trajectory,  the  search  area  was  computed  and  the  position  of  the  real  drifter  checked  to 
see  whether  it  was  inside  or  outside  the  search  area.  The  number  of  real  drifters  within  the 
search  area  was  summed  for  each  hour  and  divided  by  the  number  of  trajectories  so  that  the 
fraction  ( percent)  inside  the  predicted  search  area  could  be  computed. 

Comparison  of  the  percent  success  measure  for  the  two  turbulence  models  tested  shows  clearly 
that  the  random  flight  method  is  superior  in  both  regions  (Figure  D-7).  Approximately  80-90 
percent  of  all  drifter  search  areas  using  this  method  enclosed  the  real  drifter  position  and  there 
was  little  variation  with  time.  The  comparison,  using  the  screened  drifter  subset  is  improved 
somewhat,  peaking  above  90  percent  in  both  regions.  Search  areas  estimated  using  the  random 
walk  model  of  turbulence  are  severely  underestimated  with  percent  success,  for  the  unscreened 
subset,  dropping  from  about  70-80  percent  at  a  prediction  time  of  1  hour  to  35-60  percent  after 
24  hours.  As  for  the  random  flight  results,  the  comparison  using  the  screened  subset  improves 
typically  by  about  5  percent  at  all  prediction  times.  The  discrepancy  between  the  search  areas 
produced  using  the  two  methods  is  most  apparent  for  the  simulations  on  the  New  Jersey  shelf. 
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Figure  D-7.  Comparison  of  the  uncertainty  bounds  for  predicted  drifter  position  using  the 
random  flight  (green)  and  random  walk  (red)  turbulence  models  for  the  Block 
Island  region  (left)  and  New  Jersey  shelf  region  (right).  The  solid  curves  are 
computed  using  only  those  trajectories  that  both  start  and  end  within  the  10 
percent  coverage  zone  while  the  dashed  curves  are  computed  based  on  all 
trajectories  that  start  within  the  nominal  coverage  zone.  For  each  prediction 
time  we  plot  the  percent  of  cases  where  the  actual  drifter  location  at  that  time 
fell  within  the  estimated  95  percent  confidence  region. 


D.4  Discussion 

Drifter  positions  predicted  over  24  hours  using  CODAR  surface  currents  are  clearly  superior  to 
the  persistence  forecast  in  estimating  the  final  drifter  location.  The  mean  separation  between 
predicted  and  observed  drifter  location  is  50-70%  of  the  separation  using  the  persistence  forecast. 
Using  currents  from  a  numerical  circulation  model,  Thompson  et  al.  (2003)  simulated  the 
trajectories  of  a  number  of  surface  drifters  on  the  Scotian  shelf.  They  estimated  the  50th 
percentile  separation  value  after  24  hours  to  be  6  km,  which  is  very  similar  to  the  mean 
separation  of  6-7  km  found  in  the  present  study.  This  suggests  that  trajectory  predictions  using 
CODAR  surface  currents  have  comparable  skill  in  predicting  target  trajectories  as  predictions 
using  numerical  model  currents. 

The  search  area  evaluations  can  be  used  to  assess  the  consistency  of  the  turbulence  models  and 
their  associated  parameters.  A  consistent  turbulence  model  would  be  expected  to  provide,  for 
instance,  a  95  percent  confidence  region  that  is  successful  (that  is  including  the  real  drifter 
location)  95  percent  of  the  time.  The  random  flight  method  provides  search  areas  that  enclose 
the  real  drifter  approximately  90  percent  of  the  time,  whereas  the  random  walk  formulation  does 
significantly  worse.  This  suggests  that,  with  the  parameter(s)  estimated  from  the  drifters,  the 
random  flight  model  is  nearly  self-consistent  while  the  random  walk  model  is  definitely  not  self- 
consistent.  The  difference  between  the  random  flight  success  rate  and  the  expected  95  percent 
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suggests  a  slight  underestimate  of  either  the  variance  or  the  time  scale.  The  latter  seems  the 
more  likely  candidate  here,  as  we  have  estimated  the  turbulent  time  scale  as  half  the  value  of  the 
first  zero  crossing  of  the  autocovariance  function  and  not  the  true  integral  time  scale.  Note  that 
evaluation  of  the  integral  time  scale  by  integration  of  the  autocovariance  to  infinite  lags  is 
problematic  in  the  presence  of  the  quasi-periodic  motions  evident  in  the  covariance  functions  of 
Figure  D-3. 

The  comparison  of  the  effectiveness  of  random  walk  and  random  flight  derived  drifter  dispersion 
is  based  on  the  case  of  equivalent  turbulent  velocity  fluctuation  variance.  Alternatively,  one 
could  specify  the  diffusion  coefficient  and  either  the  fluctuation  variance  or  the  turbulent  time 
scale  for  the  random  flight  model.  In  the  case  of  equivalent  Kx,  the  results  of  Zambianchi  and 
Griffa  (1994)  show  that,  for  homogeneous,  stationary  turbulence,  the  random  walk  model 
overestimates  the  particle  dispersion  for  t  <  Tu.  At  times  large  compared  to  the  turbulent  time 
scale,  the  two  models  predict  that  the  particle  cloud  dispersion  increases  at  the  same  rate  but  the 
offset  introduced  by  the  initial  overestimate  persists.  For  turbulent  time  scales  of  1-3  hours  as 
determined  in  the  present  study,  the  difference  is  10-20  percent  at  times  after  prediction  of  12-24 
hr.  One  could  clearly  increase  the  diffusion  coefficient  in  the  random  walk  case  to  produce  more 
dispersion  of  the  pseudo-drifter  cloud,  and  if  the  diffusion  coefficient  were  set  to  the  random 
flight  value  (  crTu)  we  would  expect  that  the  search  areas  would  be  at  least  as  effective. 

D.5  Summary  and  Conclusions 

Comparison  of  real  drifter  trajectories  and  trajectories  predicted  using  CODAR-derived  surface 
currents  illustrates  the  value  of  these  data  for  search  and  rescue  operations.  For  prediction  times 
of  1-24  hr,  the  mean  (and  95th  percentile)  distance  between  the  CODAR-predicted  position  and 
the  real  position  is  smaller  than  the  distance  traveled  by  the  drifter.  This  indicates  that 
predictions  using  CODAR  velocities  are  more  accurate  than  the  so-called  “persistence”  forecast 
(zero  drifter  velocity).  Although,  not  shown  here,  CODAR  trajectory  predictions  are  also 
superior  to  those  produced  using  current  CG  practice,  in  which  the  advective  velocity  is  obtained 
from  NOAA  tidal  current  predictions  in  nearshore  waters  and  a  surface  current  climatology 
offshore  (O’Donnell,  2005). 

The  statistics  of  the  combination  of  sub-gridscale  velocity  and  CODAR  velocity  error  that 
contribute  to  the  dispersion  of  a  cloud  of  pseudo-drifters  have  been  estimated  using  the 
ensemble-averaged  covariance  functions  of  CODAR-drifter  velocity  differences.  Approximate 
consistency  of  the  estimates  of  turbulent  velocity  variance  and  timescale  was  demonstrated  for 
the  random  flight  turbulence  model  by  evaluation  of  the  resulting  search  areas,  defined  as  the 
region  in  which  95  percent  of  the  pseudo-drifters  are  located.  The  random  flight  search  areas 
include  the  real  drifter  location  in  80-90  percent  of  cases.  Using  the  turbulent  velocity  variance 
estimated  from  the  zero-lag  autocovariance  function,  as  for  the  random  flight  model,  the  random 
walk  search  areas  were  significantly  less  effective.  This  is  due  to  the  fact  that  the  effective 
turbulent  time  scale  in  the  random  walk  formulation  is  dt/2  or  0.5  hr  with  the  time  step  used  in 
this  study.  One  could  expect  to  achieve  satisfactory  search  area  predictions  using  the  random 
walk  model  by  specifying  a  diffusion  coefficient  equal  to  the  random  flight  value. 
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Abstract 

This  report  summarizes  a  hierarchy  of  stochastic  particle  models,  where  the  object’s  position, 
velocity,  and  acceleration  are  progressively  represented  as  Markovian  processes.  Numerical 
implementation  and  testing  of  the  first  and  second  order  models  (random  walk  and  random 
flight)  against  an  analytic  solution  to  the  diffusion  equation  show  very  good  agreement  provided 
that  2,000  or  more  independent  simulations  are  ensemble-averaged.  The  random  flight  model  is 
shown  to  predict  smaller  search  areas  than  the  random  walk  model,  with  a  long-term  reduction  in 
the  area  proportional  to  the  dispersion  coefficient  times  the  velocity  autocorrelation  time  scale. 
This  offset  ramps  in  immediately  after  the  release,  and  occurs  over  the  velocity  autocorrelation 
time  scale.  The  particle  models  were  applied  to  predict  the  trajectories  of  seven  CG,  Self 
Locating  Datum  Marker  Buoys  (SLDMB)  Argos  tracked  drifters  (Davis  like)  deployed  in  three 
clusters:  one  located  in  western  Block  Island  Sound,  and  the  other  two  near  the  coast  and  shelf 
break  in  the  New  York  Bight.  The  buoys  were  deployed  for  a  35-day  period  starting  on  July  27, 
2004.  High  frequency  coastal  radar  (CODAR)  measurements  were  collected  during  the  same 
time  period  by  a  short  range  system  (50  km  range  and  5  km  resolution)  operated  by  the 
Universities  of  RI  and  Connecticut  for  the  Block  Island  Sound  and  adjacent  shelf  area,  and  by 
the  long  range  system  (150  km  range,  6  km  resolution)  operated  by  Rutgers  University  for  the 
Mid  Atlantic  Bight.  The  motion  of  the  buoys  was  dominated  by  a  mean  southwesterly  shelf 
transport  and  inertial  and  semi-diurnal  tidal  oscillations.  A  cluster  analysis  of  the  very  limited 
number  of  SLDMBs  gives  dispersion  coefficients  in  the  range  of  40  to  80  m2/sec.  Analysis  of 
the  CODAR  velocity  errors  and  variances  gives  values  in  the  range  of  40  to  700  m2/sec,  with 
velocity  autocorrelation  time  scales  in  the  range  of  4  to  7  hours,  depending  on  the  velocity 
component,  the  location,  and  whether  the  current  record  used  to  determine  the  autocorrelation 
time  scale  is  de-tided.  Comparison  of  the  velocities  derived  from  the  drifters  and  the  radar 
system  shows  differences  comparable  to  the  observed  speeds.  Scatter  plots  for  the  Block  Island 
Sound  (Mid  Atlantic  Bight)  show  a  correlation  of  less  than  0.1  (greater  than  0.75)  for  the 
east/west  component  and  greater  than  0.6  (less  than  0.4)  for  the  north/south  component. 
Correlation  coefficients  were  observed  to  be  much  lower  in  areas  where  the  percent  data  return 
was  below  50  percent.  Statistically  independent  simulations  were  performed  using  SARMAP,  a 
search  and  rescue  model,  to  predict  the  daylong  trajectories  at  successive  (non-overlapping) 
locations  along  the  paths  of  the  seven  SLDMBs  in  the  Mid  Atlantic  Bight.  Model  predicted  and 
observed  locations  at  the  end  of  one  day  (typical  search  and  rescue  model  prediction  time  scale) 
showed  an  averaged  difference  of  9.3  km,  with  an  error  in  location  comparable  (about  80 
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percent)  to  the  distance  the  drifter  traveled  in  one  day.  Errors  were  significantly  higher  in  areas 
where  the  CODAR  data  return  rates  were  less  than  50  percent.  Estimates  of  the  dispersion 
coefficients  necessary  as  input  to  the  search  and  rescue  model,  to  ensure  that  the  predicted  and 
observed  locations  of  the  drifters  were  within  the  model  predicted  search  area,  were  made  and 
gave  values  ranging  from  20  to  500  m  /sec,  with  a  median  value  of  90  m  /sec.  The  largest 
dispersion  coefficients  were  associated  with  areas  where  the  radar  data  returns  were  lowest.  The 
various  estimates  of  the  dispersion  performed  (cluster  analysis,  radar  velocity  variance,  and 
trajectory  analysis)  gave  values  consistent  with  Okubo’s  relationship.  Random  flight  models 
offered  no  improvement  in  predictive  performance  for  any  of  the  applications  due  to  the 
substantial  uncertainty  in  estimates  of  the  dispersion  coefficient. 

E.l  Introduction 

In  planning  search  and  rescue  operations,  it  is  common  practice  to  use  Lagrangian  (particle) 
trajectory  models  to  predict  the  transport  of  objects  (i.e.  vessel,  person  in  the  water,  life  raft,  etc.) 
at  sea.  The  majority  of  the  models  employ  random  walk  (Markovian  in  position)  techniques  to 
perform  these  simulations  (Breivik  et  al,  2004;  Spaulding  and  Howlett,  1996;  Spaulding  and 
Jayko,  1991;  ASA,  2003).  The  current  and  wind  fields  used  to  force  these  models  are  derived 
from  hindcasts,  nowcasts,  and  forecasts,  depending  on  the  application,  and  are  supplied  by 
hydrodynamic  and  meteorological  models  or  the  direct  use  or  interpolation/extrapolation  of 
observations.  The  movement  of  the  objects  is  typically  approximated  as  the  vector  sum  of  the 
current  field  plus  an  empirically  based  down  and  cross  (leeway)  drift  in  response  to  wind  forcing. 
Allen  and  Plourde  (1999)  and  Allen  (1999)  provide  empirical  leeway  drift  factors  (drift  rate, 
down  and  cross  wind,  versus  wind  speed)  and  the  associated  uncertainties  for  63  separate  classes 
of  objects. 

The  principal  goal  of  search  and  rescue  model  simulations  is  to  predict  the  location  and  size  of 
the  area  where  the  lost  object  is  most  likely  to  be  found  (search  area)  as  a  function  of  time. 
Simulation  periods  typically  extend  from  several  days  in  the  past  to  1  or  2  days  in  the  future. 

The  search  planner  is  typically  faced  with  making  critical  assumptions  about  the  initial  release 
location  and  the  conditions  that  existed  at  that  time.  (When  and  where  did  the  vessel  experience 
its  initial  distress?  Was  the  vessel  upside  down?  Does  the  life  raft  have  a  drogue  deployed?). 
Stochastic  simulations  are  then  performed  to  predict  the  location  of  the  search  area.  Once  the 
simulation  is  complete  and  the  search  area  estimated  the  search  planner  must  determine  how  to 
allocate  search  resources  (i.e.  boats,  helicopters,  airplanes)  to  maximize  the  probability  of 
success  in  the  search.  Wind  and  current  induced  advective  processes  can  reasonably  predict  the 
center  of  the  search  area.  The  size  of  the  search  area  however,  depends  on  the  evolution  of  the 
ocean  turbulence  fields  and  the  shear  in  the  mean  current,  normally  parameterized  in  terms  of  the 
horizontal  dispersion  coefficient,  over  the  duration  of  the  simulation.  If  estimates  of  the  search 
area  are  too  small,  the  probability  of  detecting  the  object  is  reduced.  On  the  other  hand,  if  the 
search  area  is  too  large,  it  is  likely  that  search  resources  will  not  be  available  to  cover  the  search 
area,  thus  lowering  the  probability  of  detection,  since  some  of  the  area  will  remain  unsearched. 
Application  of  higher  order  stochastic  models  (e.g.  random  flight,  Markovian  in  both  position 
and  displacement)  offer  the  potential  to  reduce  the  size  of  the  search  area  because  the  effective 
dispersion,  compared  to  its  random  walk  counterpart,  is  reduced.  In  spite  of  the  potential 
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improvements,  none  of  the  current  operational  search  and  rescue  models  currently  employ  higher 
order  stochastic  trajectory  models. 

In  addition,  the  rapid  development  and  implementation  of  high  frequency,  short  and  long  range 
coastal  radar  systems  for  real-time  measurements  of  surface  currents  provide  a  new  source  of 
near-real-time  observations  of  the  surface  current  field  that  could  be  used  as  input  to  search  and 
rescue  models.  These  systems  have  the  potential  of  providing  high  spatial  (1.5  to  6  km)  and 
temporal  (1  to  3  hr)  resolution  observations  of  surface  current  fields  over  ranges  of  50  to  150  km. 
In  addition,  these  systems  can  provide  estimates  of  the  surface  velocity  variances  that  might  be 
exploited  to  improve  estimates  of  the  effective  turbulent  dispersion,  and  hence  improve 
predictions  of  the  size  of  the  search  area. 

The  first  objective  of  this  paper  is  to  summarize  a  hierarchy  of  stochastic  particle  models  from 
the  literature  that  might  be  used  in  search  and  rescue  models  to  predict  the  search  area.  The 
second  objective  is  to  apply  random  walk  and  flight  models,  driven  by  surface  current  data 
derived  from  short-  and  long-range  high  frequency  radar  systems,  to  understand  the  transport  and 
dispersion  of  seven,  7/10  scale  Davis  (1985)  like  surface  drifters  deployed  in  Block  Island  Sound 
and  the  Mid  Atlantic  Bight  by  the  CG  over  a  35-day  period  in  July  and  August,  2004.  Integral  to 
both  objectives  is  an  assessment  of  the  dispersion  characteristics  and  autocorrelation  time  of  the 
velocity  variances  of  these  coastal  waters. 

Section  2  provides  a  summary  of  the  hierarchy  of  stochastic  particle  models,  drawn  from  the 
literature  on  Lagrangian  dynamics  in  the  ocean  and  atmosphere,  which  might  be  used  in  search 
and  rescue  models.  Numerical  implementation  and  testing  of  the  models  against  an  analytic 
solution  to  the  diffusion  equation  are  also  covered  in  this  section.  Application  of  random  walk 
and  flight  models  to  the  surface  drifter  data  to  estimate  the  turbulent  dispersion  and  to  assess 
their  predictive  performance  is  provided  in  Section  3.  Anticipating  the  use  of  high  frequency 
radar  surface  current  data  to  predict  the  movement  of  surface  drifters,  Section  4  provides  a 
comparison  of  velocities  derived  from  drifter  trajectories  with  corresponding  data  collected  by 
the  high  frequency  radar.  Section  5  provides  the  associated  estimates  of  dispersion  coefficients 
based  on  the  radar  data.  Predictions  of  the  drifter  trajectories  were  then  made  using  a  state  of  the 
art,  search  and  rescue  model  (SARMAP;  ASA,  2003),  forced  by  the  observed  high  frequency 
radar  current  fields  (Section  6).  Differences  between  the  model  predictions  and  observations 
were  used  to  assess  the  dispersion  coefficient  required  in  SARMAP  to  ensure  that  the  model 
predicted  search  area  included  the  observed  location  of  the  drifters.  Study  conclusions  are 
presented  in  Section  7  and  references  in  Section  8. 

E.2  Hierarchy  of  stochastic  particle  models 

The  large  number  of  scales  participating  in  the  motion  makes  modeling  of  the  movement  of 
particles  (search  and  rescue  objects)  in  the  ocean  complicated.  For  the  search  and  rescue 
problem,  these  range  from  continental  shelf  to  molecular  scales.  Since  methods  are  not  available 
to  solve  for  all  scales  simultaneously,  it  is  common  practice  (Taylor,  1921)  to  divide  the  flow 
into  a  mean  field  representing  the  large-scale  flows,  U,  and  turbulence  representing  the  smaller 
scale  flows  (u’). 
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Following  the  work  of  Griffa  (1996),  Griffa  et  al,  (1995),  and  Berloff  and  McWilliams  (2002),  it 
is  proposed  that  the  trajectories  of  objects  drifting  at  sea  be  modeled  by  a  hierarchy  of 
progressively  more  complicated  Markovian,  stochastic  particle  models.  These  models  are  simple 
examples  of  the  generalized,  stochastic,  nonlinear,  Langevin  equations  for  arbitrary  dimensions 
(Risken,  1989).  Markovian  models  are  typically  used  to  describe  processes  whose  conditional 
probability  density  at  time  t  n  is  solely  dependent  on  values  at  an  earlier  time  t  n_i.  They  have 
been  shown  to  be  well  suited  for  modeling  the  movement  of  drifters  in  the  ocean  (Dutkiewicz  et 
al,  1994;  Griffa,  1996;  Bauer  et  al,  1998,  2002;  Mariano  et  al,  2002).  A  complete  derivation  of 
the  hierarchy  of  models  is  provided  in  appendices  to  Berloff  and  McWilliams  (2002). 

The  lowest  level  model  (Model  1 ,  alternately  referred  to  as  the  zeroth  order  Markov  model) 
assumes  that  the  particle  position,  x,  is  a  Markov  variable.  This  model  assumes  that  the  scales  of 
acceleration  and  velocity  are  infinitesimal.  This  model  can  be  shown  to  give  solutions  to  the 
advective-diffusion  equation  (Csanady,  1973)  and  is  widely  known  as  a  “random  walk”  model. 
The  model  is  simple  to  implement  and  widely  used  in  oceanographic  studies.  Model  1  is  used  in 
most  search  and  rescue  models  to  make  trajectory  predictions  (Breivik  et  al,  2004;  Spaulding  and 
Howlett,  1996;  Spaulding  and  Jayko,  1991;  ASA,  2003).  Model  2,  often  referred  to  as  a 
“random  flight”  model  (Thomson,  1986;  van  Dop  et  al,  1985),  assumes  that  both  the  position,  x, 
and  the  velocity,  u,  are  jointly  Markovian.  The  scale  of  acceleration  is  once  again  assumed 
infinitesimal.  This  model  is  well-suited  to  describe  meso-scale  particle  motions  in  the  upper 
ocean  (Griffa,  1996).  Finally,  Model  3  assumes  x,  u,  and  a  are  jointly  Markovian  (Sawford, 
(1991),  Pope,  (1994)).  This  model  is  typically  used  when  particles  display  strong  looping  or 
wave-like  motion  (Mariano  et  al,  2002).  Berloff  and  McWilliams  (2002)  extend  this  internally 
consistent  hierarchy  to  equations  for  the  hyper-acceleration  or  the  time  derivative  of  the 
acceleration.  This  model  is  most  useful  if  the  particles  move  in  vortex  or  ring-like  structures 
with  very  long  lives.  This  model  is  not  considered  here  because  it  is  well  beyond  the  scope  of 
the  current  state  of  practice,  requires  extensive  supporting  data,  and  is  not  justified  since  the 
potential  areas  of  application  do  not  have  particle  trajectories  that  display  strong  wave  or  looping 
like  movements  over  the  typical  forecast  periods. 

The  governing  equations  for  each  of  the  three  Markovian  models  are  presented  below.  Only  one 
component  is  provided,  as  the  other  is  analogous.  The  following  assumptions  have  been  made 
for  the  present  application. 

1 .  The  velocity  field  is  two-dimensional.  This  assumption  is  a  natural  result  of  focusing  on 
the  transport  of  particles  floating  at  or  very  near  the  sea  surface. 

2.  Models  describe  the  motion  of  single,  independent  particles  released  in  the  flow  field;  the 
movements  of  individual  particles  are  not  correlated  with  each  other. 

3.  The  particles  are  dynamically  passive  (i.e.  they  don’t  affect  the  fluid  dynamic  motion). 

4.  The  turbulence  velocity  is  homogenous  in  space  and  stationary  in  time.  For  the  present 
application,  it  is  assumed  that  temporal  and  spatial  variations  are  accounted  for  by  the 
temporally  and  spatially  varying  flow  fields  provided  by  the  high  frequency  surface  radar 
systems  (horizontal  spatial  resolution  on  the  order  of  1  to  6  km  and  a  temporal  resolution 
of  1  to  3  hrs)  or  hydrodynamic  models. 

5.  The  two  components  of  velocity  are  independent  (i.e.  spatial  correlations  of  the  velocity 
are  zero). 

6.  Random  forcing  is  Gaussian  (or  normal). 
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E.2.1  Model  1:  Markovian  x 


The  governing  equation  for  particle  motion,  in  incremental  form,  for  this  model  is 


dx  =  Udt  +  dx’ 

(1) 

dx’=  K 1/2  dw 

(2) 

with  initial  conditions  x(t=0)  =  0.  dx  is  the  total  particle  displacement  during  time  interval  dt,  U 
(x,  t)  is  the  mean  flow,  dx’  is  the  displacement  due  to  the  turbulent  velocity  field,  and  dw(t)  is  a 
random  increment  from  a  normal  distribution  (zero  mean  and  second  order  moment  of  2dt). 

K  is  the  turbulent  diffusion  coefficient  and  is  given 

K=  gu2  T  (3) 

2  2  2 

Where  au  is  the  velocity  variance  (cm  /sec  )  and  T  is  the  turbulence  time  scale.  In  order  for 
au2=  <(dx’/dt)  2>,  T  must  be  equal  to  dt/2. 

Griffa  (1996)  presents  the  Fokker-Planck  equation  associated  with  Model  1  and  shows  that  it 
provides  a  solution  to  the  well-known,  advective  diffusion  equation  for  the  average  concentration 
of  a  constituent,  where  K  is  the  eddy  diffusion  coefficient. 

Equation  1  shows  that  the  movement  of  a  particle  is  a  result  of  contributions  from  the  mean  flow 
and  turbulence.  The  latter  is  represented  as  a  random  stochastic  process,  uncorrelated  from  one 
time  to  the  next.  The  particle  hence  receives  a  random  pulse  or  input  due  to  actions  of  the 
turbulence  field,  but  retains  no  memory  of  the  prior  turbulent  impulse  it  experienced  at  the  earlier 
time  step.  For  this  model  T,  and  hence  dt,  are  small  compared  with  temporal  variations  of  U  and 
with  the  actual  time  t.  This  model  is  therefore  applicable  for  times  when  t »  T. 

The  dispersion,  or  mean  squared  separation  of  particles,  Si  =  <(x-  <x>)  >,  where  the  <  > 
indicates  the  mean  value,  for  this  model  is  given  by 

Si2  =  2Kt  (4) 

where  the  subscript  1  refers  to  Model  1 .  The  mean  dispersion  of  particles  is  predicted  to  grow 
linearly  with  time,  in  agreement  with  the  solution  to  the  advective  diffusion  equation  (Csanady, 
1973). 

E.2.2  Model  2:  Joint  Markovian  x  and  u 

The  incremental  equations  for  particle  motion  for  Model  2  are  given  by 


dx  =  (U+u)dt 

(5) 

du  =  -(1/T)udt  +  K*  1/2  dw* 

(6) 
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the  initial  conditions  are  x(0)  =  0  and  u(0)  =  u,  where  u  is  drawn  from  a  Gaussian  distribution 
with  a  zero  mean  and  variance  of  cu2,  dw*  is  a  random  increment  with  the  same  statistical 
characteristics  as  in  Model  1  (Equation  2),  and  T  is  the  turbulent  time  scale. 

K*  is  the  diffusion  coefficient  and  can  be  expressed  in  terms  of  the  velocity  variance  and  the 
Model  1  diffusion  coefficient,  K,  by 

K*  =  au2/T  =  K/T2  (7) 

Model  2  (Equation  5)  shows  that  perturbations  of  the  velocity  field  at  each  time  increment  are 
composed  of  two  components,  a  random  impulse  (dw*)  and  a  loss  of  velocity  (momentum) 
(-(1/T)  u  dt)  as  a  result  of  the  memory  of  the  velocity  from  the  prior  time  step.  A  particle  hence 
retains  a  memory  of  its  initial  turbulent  velocity  over  a  finite  time  of  order  T. 

T  is  typically  derived  from  the  velocity  autocorrelation  function,  expressed  as 

R(x)  =  au  "2<u(t)u(t+x)>  (8) 

where  x  is  the  time  lag.  Numerous  applications  have  shown  that  the  autocorrelation  function 
decays  exponentially  for  meso-scale  oceanographic  problems  (see  Griffa,  1996)  and  hence  R  can 
be  approximated  by 

RW-e-"  (9) 

where  T  is  the  integral  time  scale  of  the  autocorrelation  function. 

Griffa  (1996)  and  van  Dop  et  al,  (1985)  present  the  Fokker-Planck  equation  for  Model  2  and 
show  that  the  equation  includes  a  term  that  constitutes  the  memory  of  the  velocity  for  a  finite 
period  of  time. 

Model  2  can  represent  transient  processes  with  times  t  <  T;  as  well  as  those  where  t  >  T,  as  in 
Model  1 .  Model  2  fails  however  at  very  small  times  since  the  acceleration  has  an  infinitesimally 
small  scale  and  Equation  5  is  discontinuous  at  each  step. 

For  the  case  of  constant  K  and  T,  Equation  6  has  an  exact  solution  (van  Dop  et  al,  1985). 

The  particle  dispersion  for  Model  2,  S22,  for  this  case  is  given  by 

S22  =  2Kt-  2KT(l-e  _t/T)  (10) 

Equation  10  shows  that  S22  has  an  exponential  like  shape.  This  solution  also  shows  that  the 
velocity  autocorrelation  must  have  the  exponential  form  given  by  Equation  9.  Comparing  Si  and 
S2,  it  is  noted  that  Model  2  predicts  a  lower  dispersion  than  Model  1  for  t  <  T.  The  dispersion 
approaches  that  given  by  Model  1  for  large  values  of  t;  however  there  is  a  permanent  offset,  with 
a  value  of  2KT. 
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Model  3:  Joint  Markovian  x,  u,  and  a 

The  governing  equations  for  Model  3,  written  in  incremental  form,  are 
dx  =  (U+u)dt 
du  =  adt 

da  =  -(1+  a  Ta/T)dt/  Ta  -udt/(T  Ta)  +  K**1/2dw** 

where 

K**  =  (1/Ta)(aa2(l+Ta/T)) 


(11) 

(12) 

(13) 


(14) 


o  a2=  a  u2/(  T  Ta)  (15) 

the  initial  conditions  are  x(0)  =  0,  u  (0)  =  u,  a  (0)  =  a  where  the  values  for  u  and  a  are  obtained 
from  a  Gaussian  distribution,  with  zero  mean  and  variances  of  cu2  and  a2a,  respectively.  K**  is 
the  diffusion  coefficient,  a2a  is  the  acceleration  variance,  and  dw  **  is  a  random  increment  with 
the  same  statistical  characteristics  as  dw  and  dw*  in  Equations  2  and  5,  respectively.  Model  3 
has  two  time  scales,  T  and  Ta,  representing  the  velocity  and  acceleration  time  scales, 
respectively. 

Assuming  that  the  velocity  and  acceleration  autocorrelations  can  be  approximated  by  exponential 
relationships  results  in 

R(x)  =  (e  T/T  -(Ta/T)  e  T/Ta)/(1 -Ta/T)  (16) 

for  the  velocity  autocorrelation  function. 

To  understand  the  physical  behavior  of  Model  3,  as  represented  by  Equations  11  -  13,  let  us 
assume  there  is  a  complete  separation  of  the  temporal  scales  of  velocity  and  acceleration  such 
that  Ta «  T.  For  long  time  lags  (x  »  Ta)  R  (x)  reduces  to  Equation  9  (Model  2).  For  time  lags 
on  the  order  of  Ta,  R  has  a  quadratic  behavior  in  the  vicinity  of  x  =  0.  The  autocorrelation  for 
acceleration  has  an  approximate  exponential  dependence  for  short  time  lags  and  is  given  by 

Ra(x)  =  e  _T/Ia  (17) 

E.2.3  Extension  to  inhomogeneous  and  non-stationary  turbulence 

When  turbulent  parameters  in  the  above  models  are  spatially  or  temporally  varying,  the  basic 
Langevin  equations  given  above  are  not  strictly  correct.  Errors  manifest  themselves  in  the  form 
of  particle  distributions  that  are  uncharacteristically  high  in  areas  of  low  turbulence.  This 
problem  has  proven  particularly  important  in  particle  transport  problems  where  the 
temporal/spatial  variations  in  the  turbulence  levels  are  strong.  The  typical  strategy  used  to 
address  this  problem  is  the  addition  of  a  drift  correction  term  to  the  velocity  for  Model  2  and  to 
the  velocity  and  acceleration  for  Model  3  (Berloff  and  McWilliams,  2002;  Thomson,  1987; 
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Sawford  and  Yeung,  2000;  Legg  and  Raupach,  1982).  This  term  typically  involves  the  spatial 
gradients  of  the  velocity  and  acceleration  variances.  As  an  example,  the  simplest  correction  for 
Model  2,  to  account  for  inhomogeneous  and  non-stationary  turbulence,  is  to  add  the  following 
term  to  the  right-hand  side  of  Equation  6. 

do^udt 

dx 

This  term  introduces  a  mean  acceleration  acting  on  the  particles  directed  toward  the  regions  of 
higher  variance.  It  serves  to  counterbalance  the  tendency  of  particles  to  concentrate  in  areas  of 
lower  variance.  This  adjustment  ensures  that  particle  distributions  are  well  mixed. 

Data  collected  by  the  high  frequency  radar  system  in  the  Block  Island  Sound  study  area  indicate 
values  of  a  u  ranging  from  100  to  400  cm  /sec  with  a  mean  value  of  about  100  cm  /sec  (Ullman 
et  al,  2003).  For  the  present  study,  mean  variances  in  the  range  of  350  to  550  cm2/sec2  were 
observed  in  Block  Island  Sound  and  the  Mid  Atlantic  Bight.  The  larger  values  of  variance  are 
observed  at  the  outer  edge  of  the  radar  footprint  and  are  a  result  of  the  poor  geometric  resolution 
in  the  area,  the  reduced  radial  coverage  at  long  ranges,  and  lower  data  recovery  rates.  The 
velocity  variance  however,  is  typically  smoothly-varying  over  the  central  portion  of  the  footprint. 
Given  this  situation,  it  will  be  assumed,  as  a  first  approximation,  that  inhomogeneous  and  non- 
stationary  effects  can  be  ignored. 

E.2.4  Numerical  implementation 

The  procedure  used  to  solve  the  governing  equations  presented  is  provided  below.  Since  the 
method  is  the  same  for  each  model,  the  approach  for  Model  2  is  presented,  because  it  illustrates 
all  the  key  steps  and  is  less  complicated  than  Model  3. 

The  strategy  is  to  solve  Equations  5  and  6  using  a  Monte  Carlo  technique.  In  this  approach 
simulations  are  performed  for  discrete  particles  and  the  results  of  multiple  simulations  are 
ensemble-averaged  to  predict  the  advection  and  dispersion  of  the  particles  as  a  function  of  time. 

Equations  5  and  6  are  expanded  in  discrete,  incremental  form  as 

x  n+1i+1j+1  =  x  i,jn  +  (U  “i,j  +  un  ij)  At  (18) 

un+1  i+i;j+i  =  u  n  i  j  -  (1/T)  u  n  ;  ;j  At  +  (2  cu  2  /T)  1/2  wn  *  (19) 


where  the  superscript,  n,  refers  to  the  time,  and  is  stepped  in  increments  of  At.  The  subscripts  i 
and  j  refer  to  the  two-dimensional  horizontal  location  of  the  velocity.  The  location  of  the  particle 
at  n+1,  x  n+1i+1j+1,  is  given  subscript  i+1,  j+1  to  note  that  the  particle  will  be  moved  horizontally 
to  a  new  position  during  the  time  increment,  At.  U  is  the  mean  velocity  at  the  location  (i,  j)  of 
the  particle  at  time  n  and  u  is  the  turbulent  fluctuating  component  of  the  velocity.  U  can  be 
obtained  directly  from  the  output  of  the  radar  as  the  mean  value  of  the  velocity  for  each  grid.  cu  2 
is  the  variance  of  the  turbulent  velocity.  This  value  is  available  at  each  grid  point  from  the  radar 
data.  It  can  be  averaged  over  time  at  one  grid  point  or  all  grid  points  can  be  averaged  over  both 
time  and  space.  T  is  the  integral  time  scale  from  the  autocorrelation  function  (Equation  8) 
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assuming  that  the  function  is  exponential  (Equation  9).  T  is  available  either  as  an  average  value 
over  the  radar  footprint  or  for  each  observation  point.  Finally,  wn  *  is  a  random  increment  from 
a  normal  distribution  with  zero  mean  and  second  order  moment  of  2dt.  Equations  18  and  19  are 
the  simplest  approximations  to  time  stepping  of  the  particle  position.  More  sophisticated 
integration  procedures  can  be  used,  including  predictor  corrector  and  4th  order  Runge-Kutta 
methods  (Press  et  al,  1992).  These  methods  become  necessary  if  the  time  steps  of  the  simulation 
become  sufficiently  large  so  that  a  particle  is  moved  more  than  one  grid  in  one  time  step. 

In  the  computation  cycle,  one  particle  is  selected  and  Equation  1 8  is  solved  to  determine  the 
location  of  that  particle  at  time  n+1.  For  the  first  step,  xn  is  assumed  to  be  zero,  or  some  other 
prescribed  value,  and  un j  is  selected  at  random  from  a  normal  distribution  of  velocity  variance. 
Equation  19  is  then  solved  to  determine  the  value  of  the  turbulent  velocity  component  at  n  +1.  In 
the  next  time  step,  this  value  is  substituted  in  Equation  18  for  un  y.  This  sequence  of 
calculations  is  repeated  to  predict  the  position  of  the  particle  in  time.  Similar  simulations  are 
performed  for  each  additional  particle.  Finally,  ensemble  averages  are  computed  using  the 
results  of  numerous  particle  simulations  to  predict  the  mean  position  and  dispersion  of  the 
particles  at  each  time  step. 

In  the  present  application,  estimates  of  mean  currents  and  associated  estimates  of  the  velocity 
fluctuations  are  available  from  the  radar  system  at  increments  of  1  to  3  hours.  The 
autocorrelation  times  are  estimated  to  be  on  the  order  of  4  to  6  hrs.  Simulation  times  are  typically 
one  to  several  days.  The  grid  spacing  for  the  Block  Island  Sound  CODAR  system  is 
approximately  1.6  km  and  about  6  km  for  the  Mid  Atlantic  Bight  system.  For  a  mean  grid  size 
of  1.6  km  (6  km)  and  a  peak  current  of  50  cm/sec,  0.89  hrs  (3.3  hrs)  are  required  for  a  particle  to 
transit  one  grid.  To  take  maximum  advantage  of  the  temporal  variations  of  the  data,  the  time 
step  should  be  less  than  1  hour.  In  addition,  to  ensure  that  the  particles  won’t  travel  more  than 
one  grid  in  one  time  step,  the  time  step  should  be  less  than  0.989  hrs  for  the  Block  Island  Sound 
grid.  We  have  therefore  selected  a  time  step  of  15  minutes  for  all  simulations.  This  time  step 
ensures  that  particles  will  not  transit  more  than  one  grid  in  one  time  step  and  allows  us  to  take 
maximum  advantage  of  all  available  data.  This  time  step  is  also  substantially  less  than  the 
integral  time  scale  of  the  autocorrelation  function  and  hence  allows  us  to  evaluate  any 
improvement  in  performance  from  considering  the  autocorrelation  of  the  velocity  field. 

The  present  work  is  restricted  to  the  application  of  Models  1  and  2  to  predict  the  trajectories  of 
particles.  Application  of  Model  3  is  premature  and  awaits  the  results  of  the  evaluation  of  Models 
1  and  2. 

E.2.5  Model  testing 

To  verify  that  the  numerical  solutions  for  both  Models  1  and  2  have  been  implemented  correctly, 
simulations  with  fixed  parameters  for  cases  with  available  analytic  solutions  were  performed. 
Model  predictions  of  the  mean  squared  separation  of  particles,  S2  =  <(x-  <x>)2>  versus  time 
were  compared  with  the  corresponding  analytic  solutions. 

For  the  cases  studied,  the  mean  flow  was  set  to  zero  (U  =  0).  The  no  flow  case  was  selected 
because  the  only  difference  between  Models  1  and  2  is  in  their  treatment  of  the  turbulent 
transport  term.  Simulations  were  performed  for  a  period  of  5  days  with  a  time  step,  dt  =  1 5 
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minutes  (900  seconds).  The  diffusion  coefficient  (Model  1)  was  assumed  at  K  =  1  m2/sec  (or  au2 
=  2K/dt  =  22.2  cm2/sec2). 

E.2.5.1  Model  1:  Markovian  x 

Six  simulations  were  performed  increasing  the  number  of  particles  used  in  the  simulation  from 
100  to  50,000  (100,  500,  1000,  2000,  5000,  and  50,000)  and  using  the  model  parameters  above. 
Ensemble-averaged  model  predictions  of  the  mean  displacement,  Si  =  <(x-  <x>)  >  are  plotted 
versus  time  and  against  the  analytic  solution  in  Figure  E-l.  Only  simulations  for  the  500,  5000, 
and  50,000  cases  are  shown  for  clarity  in  the  figures.  For  Model- 1,  the  analytic  solution  of  the 
mean  squared  separation  of  particles  is  given  by  Equation  4,  Si2  =  2Kt.  Figures  E-2  and  E-3 
show  the  percent  error  and  cumulative  error  squared  versus  time,  respectively,  for  the 
simulations  shown  in  Figure  E-l.  The  numerical  solution  approaches  its  analytic  counterpart  as 
the  number  of  particles  used  in  the  simulation  increases.  Maximum  errors  for  the  500  particle 
case  are  about  10  percent  and  decline  to  less  than  1  percent  for  50,000  particles  (Figure  E-2). 

The  cumulative  errors  show  similar  behavior.  Errors  are  observed  to  decrease  rapidly  with  the 
initial  increase  in  number  of  particles,  and  then  more  slowly  as  the  number  of  particles  increases. 
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Figure  E-l.  Time  series  of  the  mean  squared  separation  obtained  using  Model  1.  S_xxx,  xxx 
denotes  the  number  of  particles  used. 
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Figure  E-2.  Time  series  of  percent  error  for  Model  1.  Exxx,  xxx  denotes  the  number  of 
particles  used. 


Figure  E-3.  Time  series  of  the  cumulative  squared  error  for  Model  1.  E2_xxx,  xxx  denotes  the 
number  of  particles  used. 


E-l  1 


E.2.5.2  Model  2:  Joint  Markovian  x  and  u 

Simulations,  similar  to  those  above,  were  repeated  using  Model  2,  employing  the  same 
parameters  and  assuming  an  autocorrelation  time  scale  of  1  day.  The  analytic  solution  for  the 
mean  squared  separation  of  particles  versus  time  is  given  by  Equation  10,  S22  =  2Kt  -  2KT 
(1-e  t/T).  Figure  E-4  shows  S22  versus  time  for  the  model  predicted  (for  varying  number  of 
particles)  and  analytic  solution  for  this  case.  Also  shown  for  reference  is  the  analytic  solution  for 
Model  1.  Figures  E-4  and  E-5  provide  the  percent  error  and  cumulative  error  squared  versus 
time,  respectively,  for  the  simulations  shown  in  Figure  E-3.  Model  performance,  as  a  function  of 
number  of  particles,  shows  similar  trends  as  in  Model  1,  with  the  exception  that  the  high 
frequency  variations  in  the  errors  are  dramatically  reduced.  This  is  a  direct  result  of  the  reduced 
contribution  of  the  random  component  to  turbulent  velocity  fluctuations  inherent  in  the  random 
flight  model.  Using  the  random  flight  model  clearly  shows  a  reduction  in  the  dispersion  of 
particles  compared  to  the  random  walk  model.  This  reduction  is  observed  over  a  time  scale 
equivalent  to  the  autocorrelation  time.  The  two  models  show  the  same  slope  of  the  mean 
squared  separation  distance  with  time  after  the  several  days  of  simulation.  Implementation  of  the 
random  flight  model  results  in  an  effective  reduction  (offset)  of  S22  of  2KT. 

Results  from  simulations  using  both  Models  1  and  2  indicate  that  several  thousand  particles  are 
required  to  obtain  error  levels  below  5  percent.  Error  levels  with  500  particles  are  on  the  order 
of  7.5  percent.  Predictive  performance  can  be  improved  beyond  the  5  percent  level  but  only  by 
the  use  of  a  substantially  increased  number  of  particles,  which  adversely  impacts  computational 
times. 


Figure  E-4.  Time  series  of  the  mean  squared  separation  obtained  from  Model  2.  S_xxx,  xxx 
denotes  the  number  of  particles  used.  The  analytic  solutions  for  Models  1  and  2 
are  shown  for  comparison.  The  vertical  line  represents  the  turbulent 
autocorrelation  time  scale  (1  day). 
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Figure  E-5.  Time  series  of  percent  error  for  the  Model  2.  Exxx,  xxx  denotes  the  number  of 
particles  used. 


Figure  E-6.  Time  series  of  cumulative  squared  error  for  Model  2.  E2_xxx,  xxx  denotes  the 
number  of  particles  used. 
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E.2.5.3.  Parametric  experiments 

Parametric  simulations  were  performed  varying  the  diffusion  coefficient  (K)  (1,  10,  and  20 
m2/sec  (corresponding  to  au2  =  22,  222,  and  444  cm2/sec2),  respectively  for  a  15  min.  time  step) 
and  turbulent  time  scale  (T)  (0.5,  1,  and  10  days)  to  show  that  the  model  is  capable  of  accurate 
predictions  over  the  range  of  K  and  T  values  typically  observed  in  coastal  waters. 

Time  series  predictions  of  the  mean  squared  separation,  obtained  from  Model  1  and  Model  2,  for 
K=  1,10,  and  20  m2/sec,  together  with  corresponding  analytic  solutions,  are  shown  in 
Figure  E-7.  The  turbulent  time  scale  (T)  for  Model  2  was  maintained  at  1  day.  All  simulations 
used  5,000  particles.  Numerical  model  predictions  are  seen  to  be  in  very  good  agreement  with 
the  corresponding  analytic  solutions.  The  model  performance  declines  slightly  as  the  diffusion 
coefficient  increases.  The  offset  (2KT)  between  Model  1  and  2  predictions  is  clearly  shown  and 
increases  linearly  with  the  value  of  the  diffusion  coefficient  selected. 

The  sensitivity  of  Model  2  predictions  to  variations  in  the  autocorrelation  time  scale  (T  =  0.5,  1 
and  10  days)  is  shown  in  Figure  E-8.  A  diffusion  coefficient  of  K  =  1  m2/sec  and  5,000  particles 
were  used  in  these  simulations.  Once  again,  model  predictions  are  in  very  good  agreement  with 
the  analytic  solution.  The  offset  increases  linearly  with  the  value  of  the  autocorrelation  time 
scale,  as  predicted  by  the  analytic  solution. 


2 

Figure  E-7.  Time  series  of  the  mean  squared  separation  obtained  for  K=  1,10,  and  20  m  /sec, 
using  Models  1  and  2  (T=  1  day)  and  5,000  particles.  The  analytic  solutions  for 
each  case  are  also  shown. 
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Figure  E-8.  Time  series  of  the  mean  squared  separation  obtained  for  Model  2  (T=  0.5,  1,  and  10 
days),  using  K=1  m2/sec  and  5000  particles. 


Some  interesting  insight  can  be  gained  on  the  application  of  Model  2  (random  flight),  relative  to 
Model  1  (random  walk),  for  the  search  and  rescue  problem  by  studying  the  analytic  solutions  and 
assuming  that  the  model  parameters  are  constant.  Specifically,  we  note  that: 


1 .  The  predicted  mean  location  of  a  drifting  objective  versus  time  is  independent  of  the 
model  selected,  since  the  advective  field  and  the  associated  integration  to  obtain  particle 
position  is  exactly  the  same  in  both  models. 

2.  Model  2  predicts  a  reduction  in  the  effective  search  area  (measured  by  S2),  relative  to 
Model  1 ,  beginning  immediately  after  the  initiation  of  the  simulation  and  continuing  with 
time.  The  difference  in  area  becomes  constant  and  scales  as  KT.  For  t  »T,  the  area 
ratio  (Ra)  of  Model  2  to  Model  1  can  be  approximated  by  Ra  =  1-  T/t.  The  ratio  of  search 
area  between  the  two  hence  decreases  as  time  increases.  As  an  example  for  T  =  2  days, 
Ra  is  0.6  at  t  =  5  days;  0.8  at  10  days;  and  0.1  at  20  days. 


E.3  Application  of  random  walk  and  flight  models  to  SDLMB  drifter 
trajectories 

The  CG  deployed  three  clusters  of  three  SLDMBs:  one  set  in  Block  Island  Sound  and  the  other 
two  sets  at  the  northern  end  of  the  Mid  Atlantic  Bight,  in  the  general  vicinity  of  the  New  York 
Bight.  SLDMBs  are  air  deployable,  7/10  size,  Davis  (1985)  like  surface  drifters  (70  cm  draft,  1 
m  width)  that  are  GPS  tracked  every  half  hour  by  the  Argos  satellite  system.  Table  E-l  provides 
the  drifter  serial  number,  deployment  time  (EST)  and  location,  and  the  cluster  identifier  name. 
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All  drifters  were  released  on  July  27,  2004.  Figures  E-9a  and  b  show  the  drifter  trajectories  for 
the  duration  of  the  35-day  experiment  for  the  Block  Island  Sound  and  Mid  Atlantic  Bight 
releases,  respectively. 


Table  E-l .  Summary  of  SLDMBs  deployed  on  July  27,  2004  by  CG. 


Serial  # 

Time 

Position 

Cluster 

43049 

0830 

N39-55.98  W72-40.28 

no  data 

43061 

0834 

N39-55.98  W72-40.28 

AB  B  (five  days) 

43062 

0839 

N39-55.98  W72-40.28 

ABB 

38867 

0855 

N40-17.62  W73-34.40 

no  data 

43057 

0900 

N40-17.62  W73-34.40 

AB  A 

43060 

0902 

N40-17.62  W73-34.40 

AB  A 

32779 

1009 

N41-08.31  W7 1-44.72 

Rl  A 

32773 

1011 

N41-08.31  W7 1-44.72 

Rl  A 

32776 

1014 

N41-08.31  W7 1-44.72 

Rl  A 

Drifters  43049  and  38867  ceased  operation  immediately  after  deployment  and  provided  no  useful 
data.  43061  provided  data  for  the  first  5  days  only. 

During  the  35-day  period  in  which  the  drifters  were  tracked,  measurements  of  the  near-sea 
surface  currents  were  made  by  a  three- station,  short-range  CODAR  system  covering  the  western 
end  of  Block  Island  Sound  and  the  adjacent  area  south  of  Montauk  Point  and  operated  by  the 
Universities  of  Rhode  Island  and  Connecticut  (Ullman  et  al,  2003;  Ullman  et  al,  2004;  Codiga 
and  Rear,  2004)  and  by  a  two-station,  long-range  CODAR  system  operated  by  Rutgers 
University  (Josh  Kohut)  for  the  Mid  Atlantic  Bight.  The  short-range  system  has  a  resolution  of 
about  1 .5  km  and  a  range  of  50  km  (Ullman  et  al,  2003),  while  the  long-range  system  has  a 
resolution  of  6  km  and  a  range  of  150  km  (Josh  Kohut). 

All  data  were  processed  by  the  CODAR  data  analysis  program,  and  provided  hourly  observations 
of  the  mean  surface  current  and  corresponding  variances  for  each  grid  cell.  The  quality  control 
protocol  removed  data  where  the  signal  to  noise  ratio  was  less  than  5  dB  or  data  for  which  the 
radial  lines  used  to  assemble  the  velocity  estimates  crossed  at  angles  of  30  degrees  (150  degrees) 
or  smaller  (larger).  Plots  of  the  percent  data  recovery  (above  and  below  50  percent)  over  the 
experimental  period  are  provided  in  Figure  E-9.  Recovery  rates  were  generally  good  for  the 
Block  Island  Sound  system  covering  most  of  the  radar  footprint.  Recovery  rates  for  the  Mid 
Atlantic  Bight  were  much  lower,  with  about  half  of  the  radar  footprint  giving  less  than  50 
percent  data  returns.  The  Mid  Atlantic  Bight  measurements  showed  significant  diurnal 
variability  reflecting  the  impacts  of  nighttime  atmospheric  conditions  on  the  radar  performance. 

As  a  first  step  in  the  analysis,  the  three  sets  of  drifter  data  were  used  to  estimate  the  effective 
dispersion  of  the  drifters  over  time.  The  three  sets  of  cluster  releases  are  referred  to  as: 


RI_A:  3  buoys  in  Block  Island  Sound 

AB_A:  2  buoys  near  shore  New  York  Bight 

AB  B:  2  buoys  offshore  New  York  Bight 
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Wind  time  histories  at  the  National  Oceanic  &  Atmospheric  Administration  (NOAA)/National 
Data  Buoy  Center  (NBDC)  marine  buoys  on  the  Mid  Atlantic  shelf  (44017,  44025),  in  the  New 
York  Bight  Apex  (ALSN6),  and  the  entrance  to  Buzzard’s  Bay  (BUZM3),  are  given  in 
Figure  E-10.  Winds  were  spatially  coherent  throughout  the  study  area,  with  mean  speeds  in  the 
range  of  4  to  6  m/sec. 

All  three  sets  of  SLDMBs  drifted  in  a  southwesterly  direction  (generally  shore  parallel).  The 
mean  drift  rates  (Table  E-2)  were  in  the  range  of  5  to  8  km/day.  The  speeds  were  highly  variable 
however,  with  standard  deviations  on  the  order  of  60  percent  of  the  average  values.  Maximum 
rates  were  14  to  18  km/day.  There  was  no  statistically  significant  correlation  between  the  wind 
and  drifter  speed,  even  in  the  wind-forced  band.  The  impact  of  the  predominantly  semi-diurnal 
tidal  currents  is  evident  in  all  records,  with  excursions  of  about  2  km/day  and  oriented  at  35/215 
degrees  in  the  Mid  Atlantic  Bight.  The  tidal  currents  are  substantially  stronger,  with  speeds 
above  50  cm/sec  and  oriented  in  a  north/south  direction  on  the  southern  boundary  of  Block 
Island  Sound.  Inertial  oscillations  (18.6  to  19  hour  periods)  are  prevalent  in  the  Mid  Atlantic 
Bight  drifter  paths,  as  well  as  a  lower  frequency,  large-scale  looping  motion.  This  looping 
motion  is  larger  for  cluster  AB_B  than  AB_A. 

Table  E-2.  Minimum,  average,  maximum,  and  standard  deviation  of  SLDMB  drifter  speeds  for 
each  drifter. 


SLDMB 

ID  Number 

Minimum 

(km/day) 

Average 

(km/day) 

Maximum 

(km/day) 

Standard 

Deviation 

(km/day) 

32773 

0.2 

6.3 

15.6 

3.9 

32776 

1.4 

6.9 

18.4 

4.6 

32779 

0.5 

6.3 

14.2 

4.0 

43057 

0.2 

4.9 

14.6 

3.2 

43060 

0.9 

5.1 

12.2 

2.9 

43061 

0.3 

7.8 

18.3 

4.1 

43062 

(Not  examined  due  to  short  record) 
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Figure  E-9a.  Trajectories  of  SLDMB  drifters  released  in  the  Block  Island  Sound  area.  The 
drifter  numbers  are  noted.  The  contours  show  the  greater  than  50  percent  (light 
gray)  and  less  than  50  percent  (dark  gray)  data  returns  for  the  Universities  of 
Rhode  Island  and  Connecticut  CODAR  systems. 


E-18 


Figure  E-9b.  Trajectories  of  SLDMB  drifters  released  in  the  Mid  Atlantic  Bight  area.  The 

drifter  numbers  are  noted.  The  contours  show  the  greater  than  50  percent  (light 
gray)  and  less  than  50  percent  (dark  gray)  data  returns  for  the  Rutgers  University 
long  range  CODAR  system. 


Power  spectra  for  the  wind  speed  (buoy  44025)  (upper  panel)  and  the  east/west  and  north/south 
components  of  the  currents,  based  on  data  from  SLDMB  43057,  at  CODAR  cell  349  (location: 
73.88669  Longitude,  39.33987  Latitude,  cell  with  highest  data  return  in  MAB),  and  at  the 
CODAR  cells  along  the  trajectory  of  43057  (lower  panels),  are  shown  in  Figure  E-l  la.  A  more 
detailed  version  of  the  power  spectra,  noting  the  principal  diurnal,  semi-diurnal  and  over-tide 
components  is  provided  in  Figure  E-l  lb  for  buoy  43057  and  its  companion  in  cluster  AB  A, 
43061.  The  buoy  data  show  that  the  most  energy  is  found  at  the  inertial  oscillation  period, 
followed  by  that  in  the  semi-diurnal  tidal  band.  Energies  in  the  tidal  harmonics  and  at  periods 
typical  of  wind  forcing  (several  days)  (Figure  E-10)  are  substantially  smaller.  The  energy  at  the 
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inertial  oscillation  period,  based  on  the  SLDMB  and  the  corresponding  CODAR  data  at  the  same 
locations,  shows  that  the  CODAR-based  estimates  have  lower  energies  in  the  north/south 
direction  than  in  the  east/west  direction,  while  they  are  comparable  in  the  SLDMB  observations. 
The  CODAR  measurements,  hence,  under  represent  the  energy  at  the  inertial  oscillation  period 
in  the  north/south  direction.  The  energies  at  the  semidiurnal  tidal  band,  however,  are 
comparable. 
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Figure  E-10.  Wind  time  series  for  July  27  to  September  1,  2004  at  NOAA/NBDC  stations 
44017,  BUZM3,  44025,  and  ALSN6. 
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Figure  E-l  la.  Power  spectra  of  the  wind  (buoy  44025)  and  east/west  and  north/south 

components  of  the  currents  derived  from  the  SLDMB  43057,  from  CODAR 
grid  cell  349,  and  from  CODAR  along  the  trajectory  of  SLDMB  43057. 
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Figure  E-l  lb.  Power  spectra  for  the  east/west  and  north/south  components  of  the  currents 
derived  from  the  SLDMB  43057  and  43061.  The  principal  diurnal,  semi¬ 
diurnal  and  harmonic  tidal  frequencies  are  noted. 


2 

For  each  cluster,  the  mean  squared  separation  distances,  S  ,  were  estimated  by 

S2  (m2)  =  (Drifter  location  -  Center  of  the  cluster) 2  /  (Number  of  drifters) 

Figure  E-l 2  shows  the  mean  squared  separation  distance  (log  scale)  versus  time  for  each  cluster 
of  drifters.  Drifters  released  at  both  AB  A  and  AB  B  stayed  remarkably  close  together  for  the 
first  15  days,  with  separation  variances  of  less  than  1.0E6  to  5.0E6  m2  (distances  of  1  to  2.3  km). 
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Cluster  AB  A  started  to  disperse  significantly  after  25  days.  The  mean  squared  separation 
distance  increased  progressively  throughout  the  experiment  for  cluster  RI_A. 

Assuming  that  the  mean  squared  separation  was  approximated  by  either  a  random  walk 
(Equation  4)  or  a  random  flight  model  (Equation  10),  estimates  of  the  dispersion  coefficients 
were  made  for  each  cluster. 


Time  (days) 

Figure  E-12.  Mean  squared  separation  distance  squared  versus  time  for  each  cluster  of  drifters. 


Figures  E-13a  and  b  show  the  mean  squared  separation  distance  versus  time  for  cluster  RI-A  and 
AB  B,  respectively.  Results  for  cluster  AB_B  are  not  provided,  since  one  drifter  only  provided 
data  for  the  first  5  days  of  the  experiment.  Lower,  mid,  and  upper  bound  fits  of  the  random  walk 
model  to  the  data  were  made  and  resulted  in  the  dispersion  values  shown  in  Figure  E-13  and 
Table  E-3.  Restricting  attention  to  the  median  values,  cluster  RI_A  gave  a  dispersion  coefficient 
of  40  m2/sec,  consistent  with  what  one  might  expect  on  the  shelf.  The  values  for  AB  A  (85 
m2/sec)  were  about  twice  those  for  RI_A.  A  closer  look  at  the  mean  squared  separation  distance 
for  AB_A,  however,  shows  that  the  dispersion  estimate  is  dominated  by  the  last  10  days  of  the 
record.  If  attention  is  restricted  to  the  initial  five  days  for  AB  A  and  AB  B  (Figure  E-14),  the 
dispersion  coefficient  is  very  low,  on  the  order  of  0.3  m2/sec. 
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Table  E-3.  Dispersion  values  estimated  from  cluster  separation  statistics  assuming  a  random 
walk  model. 


Cluster 

Dispersion  coefficient  (m2/sec) 

RI_A 

AB_A 

AB_A &  AB_B* 

Lower: 

17 

20 

0.1 

Median: 

40 

85 

0.3 

Upper  bound: 

110 

185 

0.75 

For  comparison  to  the  values  in  Table  E-3,  Okubo’s  (1971)  4/3rds  power  law  gives  dispersion 
coefficients  in  the  range  of  20  to  200  m2/sec  for  time  scales  from  5  to  30  days.  These  values  are 
broadly  consistent  with,  but  larger  than,  those  observed  by  the  cluster  analysis.  Elliott  et  al, 
(1997)  have  recently  analyzed  data  from  short-term  (less  than  2  days)  dye  dispersion  studies  in 
coastal  and  estuarine  waters,  and  find  that  the  dispersion  coefficient  scales  linearly  with  the 
current  speed  for  these  releases.  For  peak  tidal  current  speeds  in  the  range  of  40  to  60  cm/sec  for 
the  present  study  areas,  dispersion  coefficients  would  be  approximately  0.4  to  0.6  m2/sec. 

It  is  impossible  to  reconcile  the  differences  in  estimates  of  the  dispersion  coefficient  from  the 
cluster  analysis  given  the  very  limited  drifter  data  available. 


Time  (days) 

Figure  E-13a.  Mean  squared  separation  distance  versus  time  for  cluster  RI_A.  Also  shown  are 
lower,  mid  and  upper  bound  fits  of  random  walk  and  random  flight  models  with 
the  parameter  values  given  in  the  legend. 
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Figure  E-13b.  Mean  squared  separation  distance  versus  time  for  cluster  AB_A.  Also  shown  are 
lower,  mid,  and  upper  bound  fits  of  random  walk  and  random  flight  models,  with 
the  parameter  values  given  in  the  legend. 
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Figure  E-14.  Mean  squared  separation  distance  versus  time  for  cluster  AB  A  and  AB-B 
restricted  to  the  first  5  days  of  the  experiment.  Also  shown  are  lower,  mid, 
and  upper  bound  fits  of  a  random  walk  model  to  the  data,  with  the  values  of 
the  dispersion  coefficient  given  in  the  legend. 


In  order  to  apply  the  random  flight  model  to  the  data,  an  estimate  of  the  autocorrelation  time 
scale  is  required.  The  autocorrelation  function  versus  lag  time  for  the  east/west  and  north/south 
components  of  the  velocity  was  estimated  at  that  location  with  the  largest  number  and  best 
quality  data  return  in  the  Mid  Atlantic  Bight  footprint  (Station  349,  73.88669  Longitude 
39.33987  Latitude).  The  autocorrelation  for  the  two  components  are  shown  in  Figures  E-15a  and 
b.  The  analysis  was  performed  using  the  hourly  sampled  CODAR  data  and  with  a  corresponding 
de-tided  data  set  where  the  eight  largest  tidal  harmonics  were  removed.  The  tidal  currents  were 
dominated  by  the  M2,  S2,  and  N2  constituents,  with  mean  speeds  of  1 1 .3,  2.66  and  2.3  cm/sec  and 
directions  of  35,  57,  and  28.8  degrees,  respectively.  Assuming  that  the  autocorrelation  function 
can  be  represented  by  an  exponential  relationship  (Equation  9),  the  autocorrelation  time  scale 
(1/e  folding  time)  is  3.5  hrs  for  the  original  record  and  6.5  hrs  for  the  de-tided  record.  Even 
when  the  de-tided  record  is  used,  the  autocorrelation  function  shows  evidence  of  a  periodic 
variation.  Estimates  were  made  of  the  dispersion  coefficient  using  Equation  3,  with  T  based  on 
the  autocorrelation  time  scale.  The  east  and  north  components  gave  dispersion  coefficients  of 
170  and  106  m2/sec,  respectively,  for  the  original  record  and  187  and  130  m2/sec  for  the  de-tided 
record.  These  values  are  comparable  to  those  based  on  the  analysis  of  drifter  cluster  AB  A 
above. 

Estimates  of  the  dispersion  coefficient  can  also  be  made  by  directly  integrating  the 
autocorrelation  function  (Batchelor,  1953).  In  the  present  case,  the  estimate  proved  to  be  very 
sensitive  to  the  length  of  the  integration  time  since  the  autocorrelation  function  is  quite  variable 
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at  long  lag  times  (Figure  E-15).  This  is  particularly  true  for  the  original  record,  which  is  strongly 
impacted  by  tidal  variations.  Ojo  and  Bonner  (2000)  had  a  similar  experience  in  deriving 
dispersion  coefficients  from  HF  radar  data  collected  in  Corpus  Christi  and  Matagorda  Bay, 

Texas.  The  mean  value  of  their  dispersion  coefficients  was  in  the  range  of  1.5  to  4.5  m2/sec, 
while  peak  values  reached  several  thousand. 
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Figure  E-15a.  Autocorrelation  of  the  east  component  of  the  velocity  versus  lag  time  for  the 
original  CODAR  and  de-tided  CODAR  data  from  location  349  in  the  Mid 
Atlantic  Bight.  Also  shown  are  estimates  of  the  dispersion  coefficient  based  on 
the  assumption  that  Equation  10  can  represent  the  autocorrelation  function. 
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Figure  E-15b.  Autocorrelation  of  the  north  component  of  the  velocity  versus  lag  time  for  the 
original  CODAR  and  de-tided  CODAR  data  from  location  349  in  the  Mid 
Atlantic  Bight.  Also  shown  are  estimates  of  the  dispersion  coefficient  based  on 
the  assumption  that  Equation  10  can  represent  the  autocorrelation  function. 


With  an  estimate  of  the  autocorrelation  time  scale,  one  can  assess  the  impact  of  assuming  that  a 
random  flight  model  approximates  the  dispersion  of  the  drifters.  Restricting  attention  to  cluster 
RI_A,  Figure  E-13a  shows  a  comparison  of  a  random  flight  model  fit  to  the  mean  squared 
separation  distance  data,  assuming  a  median  dispersion  coefficient  of  40  m2/sec  and 
autocorrelation  time  scales  of  6,  12,  24,  and  48  hours.  Six  hours  is  approximately  the  value 
determined  from  an  analysis  of  the  CODAR  data  (Figures  E-15a  and  b).  It  is  noted  that  the 
variations  of  the  autocorrelation  time  scale  over  the  range  selected  have  very  little  impact  on 
estimates  of  the  mean  squared  separation  distance.  It  is  clear  that  there  is  insufficient  data  to 
estimate  the  improvements  resulting  from  use  of  a  random  flight  model.  Application  of  an 
unconstrained  two-dimensional  optimization  method  gives  a  K=102  m2/s  and  T=7.6  hr  for 
cluster  AB  A.  Once  again  the  use  of  a  random  flight  model  provides  no  benefit  over  a  random 
walk  approximation  in  the  present  application. 

E.4  Comparison  of  velocities  derived  from  SLDMB  versus  CODAR 

Anticipating  the  use  of  the  CODAR  data  as  input  to  a  search  and  rescue  model,  it  would  be 
useful  to  compare  the  observed  surface  current  data  from  the  CODAR  system  and  similar  data 
derived  from  the  SLDMB  drifters.  The  velocities  of  the  seven  SLDMB  drifters  were  therefore 
estimated  at  1/2  hour  time  intervals  by  dividing  the  distance  traveled  by  the  elapsed  time  in  the 
east/west  and  north/south  directions.  SLDMB  data  were  available  every  half  hour.  These  values 
were  extremely  erratic  due  to  ubiquitous  irregularities  in  the  original  position  versus  time  files. 
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The  data  were  edited  to  remove  outliers,  the  majority  of  which  were  not  contiguous  in  time.  The 
effect  of  the  editing  was  a  smoother  and  more  realistic  velocity  field.  It  was  assumed  that  the 
drifters  and  CODAR  data  both  represented  the  currents  in  the  upper  meter  of  the  water  column, 
and  hence  the  two  measurements  could  be  compared  directly. 

The  CODAR  velocity  at  a  particular  SLDMB  location  and  at  a  given  time  was  estimated  using 
inverse-distance  weighting  of  the  nearest  five  CODAR  grid  points  in  space  and  linear 
interpolation  in  time.  If  the  distance  to  any  CODAR  grid  point  was  greater  than  4  km,  the  data 
were  not  used. 

Figures  E-9a  and  b  show  the  SLDMB  trajectories,  as  well  as  the  CODAR  data  coverage  (50 
percent>  and  <50  percent)  contours  over  the  experimental  period  for  the  Block  Island  Sound  and 
Mid  Atlantic  Bight  areas,  respectively.  Comparison  of  the  above  and  below  50  percent  coverage 
in  the  two  areas  shows  that  the  difference  is  minimal  in  the  Block  Island  Sound  and  substantial  in 
the  Mid  Atlantic  Bight. 

Time  series  plots  were  prepared  to  compare  the  CODAR-observed  and  SLDMB-derived 
currents,  and  are  provided  in  Figures  E-16a,  b,  and  c  for  the  Block  Island  Sound  (a,  b)  and  Mid 
Atlantic  Bight  (c)  drifter  clusters,  respectively.  The  upper  panels  show  the  SLDMB  (thick)  and 
CODAR  (thin)  velocities,  and  the  lower  panel  the  corresponding  differences.  Figure  E-17 
provides  scatter  plots  for  the  east  (light  gray)  and  north  (black)  velocity  and  the  associated 
correlation  coefficients. 

In  general,  the  difference  in  velocity  between  the  drifter-derived  currents  and  those  from 
CODAR  are  similar  in  magnitude  to  their  speed.  The  variance  and  the  autocorrelation  time  of 
the  difference  (errors)  are  also  similar  for  both  SLDMB  and  CODAR  velocities.  Tidal  variations 
are  dominant  in  Block  Island  Sound  and  are  clearly  evident  in  both  data  sets.  They  are  less 
prevalent  in  the  Mid  Atlantic  Bight  observations.  For  the  three  buoys  released  in  Block  Island 
Sound  (32773,  32776,  and  32779),  the  scatter  plots  (Figure  E-17)  show  a  correlation  of  less  than 
0. 1  for  the  east  component  and  greater  than  0.6  for  the  north  component.  This  is  consistent  with 
the  dominant  axis  of  the  tidal  current  ellipses,  which  are  generally  north/south  in  the  area  of  the 
releases.  For  the  Mid  Atlantic  Bight  buoys,  the  scatter  plots  show  a  correlation  coefficient  of 
greater  than  0.75  for  the  east  component  (onshore-offshore)  and  less  than  about  0.4  for  the  north 
component.  This  behavior  is  just  the  opposite  of  the  releases  in  Block  Island  Sound.  The  radar 
obtains  a  more  accurate  measure  of  the  onshore-offshore  flow  than  the  along-shelf  flows.  Both 
correlation  coefficients  are  observed  to  be  lower  in  areas  where  the  percent  data  return  is  below 
50  percent  (compare  correlation  coefficients  for  43061  with  those  for  43057  and  43060). 
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Figure  E-16a.  Time  series  of  the  SLDMB  and  CODAR  derived  currents  (upper  three)  and  the 
difference  between  the  two  (lower  three)  for  drifters  released  in  Block  Island 
Sound  (see  Figure  E-9a). 
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Figure  E-16b.  Time  series  of  the  SLDMB  and  CODAR  derived  currents  (upper  three)  and  the 
difference  between  the  two  (lower  three)  for  drifters  released  in  Block  Island 
Sound  and  transported  into  the  Mid  Atlantic  Bight  (see  Figures  E-9a  and  b). 
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Figure  E-16c.  Time  series  of  the  SLDMB  and  CODAR  derived  currents  (upper  four)  and  the 
difference  between  the  two  (lower  four)  for  drifters  released  in  the  Mid  Atlantic 
Bight  (See  Figure  E-9b). 
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Figure  E-17.  Scatter  plots  for  the  east  (u,  light  gray)  and  north  (v,  black)  components  of  the 
velocity  derived  from  the  CODAR  and  SLDMB  data,  with  associated 
correlation  coefficient.  The  drifter  number  is  provided  in  the  top  left  comer  of 
each  plot. 
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E.5  Estimates  of  dispersion  coefficients  based  on  CODAR  data 

Estimates  of  the  dispersion  coefficients  for  the  Block  Island  Sound  and  Mid  Atlantic  Bight 
region  were  made  based  on  an  analysis  of  all  the  data  collected  over  the  SLDMB  deployment 
period,  where  the  CODAR  data  return  was  above  50  percent.  Table  E-4  summarizes  the  bulk- 
averaged,  east/west  and  north/south  velocity  variance,  average  error  variance,  and  the  median 
autocorrelation  time  scale.  The  velocity  variance  and  estimates  of  the  velocity  error  were 
obtained  directly  from  the  CODAR  data  processing  system.  The  dispersion  coefficients  were 
calculated  using  Equation  3,  assuming  that  T  is  the  autocorrelation  time  scale,  which  decreased 
exponentially  with  time.  Dispersion  coefficients  based  on  both  the  velocity  and  errors  variance 
are  provided.  Estimates  based  on  the  separation  analysis  presented  above  in  Table  E-3  are 
included  for  comparison. 

Table  E-4.  Estimates  of  the  dispersion  coefficients  based  on  CODAR  data  for  Block  Island 


Sound  and  Mid  Atlantic  Bight. 


Parameter 

Block  Island  Sound 

Mid  Atlantic  Bight 

East/west 

North/south 

East/west 

North/south 

Average  variance  (cnrf/sec") 

552.6 

438.3 

331.6 

345.9 

Average  error  (cnf /sec") 

45.3 

29.8 

107.3 

217.0 

Median  T*  (hrs) 

3.7 

2.8 

5.0 

5.0 

Dispersion  coefficient  (nf/sec) 

Based  on  variance 

726.2 

441.8 

596.8 

622.6 

Based  on  error 

59.5 

30.1 

193.2 

390.6 

Based  on  separation  analysis 

40 

85 

*  T  -  autocorrelation  time  scale 

Using  the  velocity  variances  gives  very  large  values  of  the  dispersion  coefficient  for  both  areas, 
particularly  for  Block  Island  Sound.  These  large  values  are  at  least  partially  attributable  to  the 
horizontal  shear  in  the  tidal  current  fields  in  Block  Island  Sound.  Replacing  the  variance  with 
the  CODAR  velocity  errors  and  using  the  autocorrelation  time  of  the  velocity  reduces  the 
dispersion  to  the  range  of  separation  analysis-based  values. 

E.6  SARMAP  simulations  of  SLDMB  trajectories 

SARMAP  (ASA,  2003)  simulations  were  performed  using  the  CODAR  data  as  input  to  predict 
the  movement  of  the  trajectories  of  the  SLDMB  drifters  for  the  Mid  Atlantic  Bight  releases  (and 
for  the  Block  Island  releases  as  they  drifted  into  the  area  covered  by  the  Mid  Atlantic  Bight  radar 
system).  Each  day-long  simulation  started  at  the  reported  SLDMB  location  at  midnight.  A  1- 
day  simulation  interval  was  used,  as  it  is  typical  of  the  time  scale  of  a  search  forecast  period. 
Simulations  were  repeated  with  a  new  start  at  the  current  SLDMB  position  each  succeeding  day, 
with  no  overlap  in  time.  Based  on  the  availability  of  SLDMB  data,  178  independent  simulations 
(Table  E-5)  were  performed. 
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Table  E-5  .  Number  of  simulations  performed  for  each  SLDMB. 


SLDMB 

Period  (Y2004) 

Simulations 

32773 

Aug  2  -  Aug  31 

30 

32776 

Aug  2  -  Aug  31 

30 

32779 

Aug  2  -  Aug  31 

30 

43057 

Augl  -  Aug  31 

31 

43060 

Augl  -  Aug  31 

31 

43061 

Augl  -  Aug  31 

26 

43062 

(not  used  as  the  record  is  too  short) 

Total 

178 

This  series  of  simulations  was  performed  to  characterize  the  advection  transport  fields  and  hence, 
no  diffusion/dispersion  was  used.  A  typical  simulation  ended  with  the  predicted  location  of  the 
drifter  some  distance  away  from  the  observed  location.  The  difference  in  final  location  between 
the  observed  and  model  predictions  is  a  result  of  errors  in  the  buoy  location  and  the  CODAR 
estimated  currents  and  the  turbulence  inherent  both  in  the  buoy  movements  and  CODAR 
observations.  Independent  of  the  reasons,  the  difference  in  location  provides  an  estimate  of  the 
equivalent  diffusion  value  such  that  the  observed  location  of  the  buoy  was  within  the  model 
predicted  search  area.  With  observed  velocities  in  range  of  5  to  18  km/day  and  CODAR  grid 
resolutions  of  6  km  in  the  Mid  Atlantic  Bight,  trajectory  predictions  sampled  no  more  than  a  few 
grids  of  current  data  in  the  day-long  simulations;  hence,  current  shear  cannot  be  an  important 
contributing  factor  to  the  model-predicted  position  of  the  drifters. 

Figure  E-18  shows  the  trajectories  of  each  buoy  (continuous  dark  gray  thick  line)  and  the 
corresponding  model  predicted  trajectories  (black  lines),  at  1-day  intervals.  A  circle  is  drawn  for 
each  simulation  centered  at  the  model-predicted  buoy  location  at  the  end  of  one  day  with  a  radius 
equivalent  to  the  distance  from  the  predicted  to  the  observed  location  at  the  end  of  the  day.  This 
circle  represents  the  minimum  area  required  so  that  the  target  (search  object)  is  within  the  model- 
predicted  search  area.  The  above  and  below  50  percent  CODAR  data  return  are  noted  by  the 
contour  intervals.  In  general,  the  circle  sizes  approximately  inversely  reflect  the  CODAR 
percent  coverage,  the  higher  the  coverage  the  smaller  the  circle  size.  The  circle  sizes  are 
particularly  large  immediately  south  of  Long  Island.  When  CODAR  returns  are  higher,  better 
estimates  of  the  velocity  are  obtained  and  reflected  in  lower  differences  between  predictions  and 
observations. 
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Figure  E-18.  Observed  SLDMB  trajectories  and  corresponding  SARMAP  predictions  for  1-day 
simulations,  restarted  every  day  along  the  SLDMB  path.  Also  shown  is  an 
inclusion  circle,  which  is  bounded  by  the  ends  of  the  observed  and  model 
predicted  trajectories. 


Figure  E-19  provides  histograms  of  the  separation  distance  between  the  predicted  and  observed 
locations  of  the  buoys  at  3,  6,  12,  and  24  hrs  after  the  start  of  the  simulation.  The  last  panel  is  a 
repeat  of  the  results  at  24  hrs  showing  additional  details  of  the  frequency  distribution  as  well  as 
the  cumulative  distribution.  Figure  E-20  shows  the  same  data  but  in  terms  of  the  average  and 
95th  percentile  separation  distances  versus  time.  Both  plots  show  that  the  separation  distance 
increases  and  the  distribution  gets  broader  and  flatter  with  time.  At  the  end  of  the  simulations, 
the  maximum  separation  distance  is  ~30  km.  The  maximum  occurrence  distance  is  in  the  range 
of  3  to  4  km.  The  median  distance  is  7  km  and  the  average  is  9.3  km.  The  95th  percentile 
separation  distance  is  23  km.  These  results  are  consistent  with  similar  analyses  performed  on 
CODAR  data  and  SDLMB  buoy  deployment  in  Middle  Atlantic  Bight  by  Ullman  et  al,  (2002), 
which  show  95th  percentile  separation  distances  at  the  end  of  1  day  of  25  km. 
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Figure  E-19.  Histograms  of  the  separation  distance  between  model-predicted  and  observed  for 
the  SLDMBs  at  3,  6,  12,  and  24  hrs  after  the  start  of  the  simulation.  The  last 
panel  shows  the  distribution  and  cumulative  separation  distance  distribution  at  the 
end  of  1  day. 
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Figure  E-20.  Average  and  95th  percentile  separation  distance  as  a  function  of  time  from  the  start 
of  the  simulation. 

Figure  E-21  shows  a  plot  of  the  separation  distance  at  the  end  of  the  day  versus  the  SLDMB  drift 
speed  for  that  day.  In  general,  the  trend,  while  statistically  very  weak,  shows  that  the  separation 
distance  increases  as  the  mean  drift  speed  increases.  The  separation  distance  is  approximately 
80%  of  the  observed  buoy  movement  over  1  day. 
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Figure  E-21.  Separation  distance  between  predicted  and  observed  locations  of  the  SLDMB  at 
the  end  of  each  day  versus  distance  traveled  by  the  SLDMB  during  that  day.  A 
linear  curve  fit  to  the  data  is  provided  (solid  line)  with  the  corresponding 
correlation  coefficient. 


Using  Equation  4  (random  walk  model)  and  assuming  that  t  is  the  time  from  release  and  the 
separation  distance  is  equivalent  to  Si,  an  estimate  can  be  made  of  the  dispersion  coefficient  that 
must  be  used  in  the  search  and  rescue  model  in  order  to  ensure  that  model  predictions  give  a 
search  area  that  includes  the  observed  location  of  the  buoy.  Figure  E-22  shows  a  histogram  of 
the  dispersion  coefficients  based  on  this  analysis. 
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Figure  E-22.  Histogram  of  dispersion  coefficients  based  on  an  application  of  a  random  walk 
model  to  the  separation  distances  at  the  end  of  1-day  simulation. 


A  diffusion  coefficient  of  40  m  /sec  is  the  most  frequent  match  (-20  percent)  and  represents  up 
to  35  percent  of  the  simulations  and  100  m2/sec  up  to  60  percent.  The  larger  values  of  dispersion 
(greater  than  500  m2/sec)  occur  in  areas  of  low  radar  data  returns  and  are  not  representative.  The 
values  of  the  dispersion  coefficients  from  this  estimate  are  consistent  with  those  estimated  using 
the  SLDMB  cluster  analysis  above. 

If  the  diffusion  coefficient  were  kept  the  same  and  the  random  flight  were  employed,  the  search 
radius  would  be  expected  to  shrink  and  miss  the  target  depending  on  the  autocorrelation  time,  or 
conversely,  the  diffusion  coefficient  must  be  larger  than  the  values  shown  in  Figure  E-22. 

E.7  Conclusions 

A  review  of  stochastic  particle  models  is  presented  and  shows  a  well-established  hierarchy  with 
position,  velocity  and  acceleration  progressively  described  as  Markovian  processes.  The 
governing  equations  for  each  are  provided.  Numerical  implementation  and  testing  against  an 
analytic  solution  for  diffusive  spreading  of  both  random  walk  and  flight  models  are  presented. 
The  numerical  model  is  shown  to  give  predictions  in  good  agreement  with  the  analytic  solution 
when  the  results  are  based  on  2,000  or  more  independent  simulations.  Parameteric  studies  show 
that  the  random  flight  model  predicts  a  smaller  search  area,  compared  to  the  random  walk  model 
beginning  immediately  after  the  initiation  of  the  simulation  and  continuing  with  time.  The 
difference  in  area  becomes  constant  and  scales  as  the  diffusion  coefficient,  K,  times  the  velocity 
autocorrelation  time  scale,  T.  For  time  »T,  the  search  area  ratio  (Ra)  can  be  approximated  by 
Ra  =  1-  T/t.  The  ratio  of  the  search  areas  hence  decreases  as  time  increases.  The  random  flight 
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model  is  therefore  expected  to  have  the  largest  impact  in  reducing  the  size  of  the  search  area  at 
time  scales  shorter  than  T. 

Based  on  the  drifter  cluster  analysis,  dispersion  coefficients  of  40  to  85  m  /sec  were  estimated 
for  the  35-day  experiments  by  fitting  the  observations  to  a  random  walk  model.  The  variability 
in  these  estimates  was  very  large  however,  covering  a  range  from  20  to  almost  200  m2/sec. 

These  dispersion  coefficients  are  broadly  consistent  with  Okubo’s  (1971)  4/3rds  power  law  and 
within  the  range  of  variability  of  his  curve  fit  to  field  data.  The  high  level  of  uncertainty  in  the 
estimates  is  clearly  a  result  of  the  very  limited  number  of  drifter  trajectories  available  to  perform 
the  analysis.  One  surprising  finding  was  that  several  of  the  drifters  in  the  Mid  Atlantic  Bight 
region  tracked  very  closely  with  one  another  over  a  period  of  almost  25  days.  The  upper  bound 
separation  distances  were  in  the  range  of  1  to  3  km.  This  implies  very  low  dispersion  (0.5 
m2/sec)  for  this  particular  cluster  of  drifters. 

Random  flight  models  were  applied  to  the  drifter  data  to  determine  their  impact  on  dispersion  of 
the  drifters.  The  autocorrelation  time  scale  was  estimated  using  CODAR  data  and  gave  values  in 
the  range  of  3.5  to  6.5  hrs  for  the  Mid  Atlantic  Bight  when  the  original  and  de-tided  data, 
respectively,  were  used.  No  improvement  was  observed  in  the  ability  of  a  random  flight  model 
to  predict  the  dispersion  of  the  drifters  because  of  the  large  uncertainty  in  the  dispersion 
coefficient.  Based  on  the  present  results,  there  is  no  motivation  to  implement  higher  order 
stochastic  trajectory  models  for  search  and  rescue  until  better  estimates  of  turbulent  dispersion 
become  available. 

Direct  comparisons  of  the  velocities  derived  from  the  drifters  and  CODAR  observations  showed 
that  the  differences  are  similar  in  magnitude  to  their  speed.  The  variance  and  the  autocorrelation 
time  of  the  difference  (errors)  are  also  similar.  Semi-diurnal  tidal  variations  are  clearly  evident 
in  both  data  sets;  dominant  in  Block  Island  Sound,  but  much  smaller  in  the  Mid  Atlantic  Bight. 
Inertial  oscillations,  with  periods  of  about  19  hrs,  are  dominant  in  the  Mid  Atlantic  Bight,  but  not 
observed  in  Block  Island  Sound.  Scatter  plots  for  the  Block  Island  Sound  cluster  show  a 
correlation  of  less  than  0.1  for  the  east/west  component  and  greater  than  0.6  for  the  north/south 
component.  For  the  Mid  Atlantic  Bight  clusters  the  scatter  plots  show  a  correlation  coefficient  of 
greater  than  0.75  for  the  east/west  component  (onshore-offshore)  and  less  than  about  0.4  for  the 
north/south  component.  This  behavior  is  just  the  opposite  of  the  Block  Island  Sound  cluster. 
Correlation  coefficients  are  observed  to  be  much  lower  in  areas  where  the  percent  data  return  is 
below  50  percent. 

Estimates  of  the  dispersion  coefficients  for  both  north/south  and  east/west  components  were 
made  using  the  CODAR  variance  and  error  data  and  the  median  autocorrelation  time  scale  for 
Block  Island  Sound  and  the  Mid  Atlantic  Bight.  The  dispersion  coefficients,  based  on  the 
velocity  variances,  were  comparable  in  the  two  directions  in  both  areas  with  values  in  the  range 
of  450  to  700  m2/sec.  Using  velocity  errors,  the  dispersion  coefficients  were  in  the  range  of  30  to 
60  m2/sec  in  Block  Island  Sound  and  200  to  400  m2/sec  in  the  Mid  Atlantic  Bight.  These 
compare  to  median  values  of  40  and  85  m2/sec  for  Block  Island  Sound  and  Mid  Atlantic  Bight, 
respectively,  based  on  the  cluster  analysis. 

A  total  of  178  statistically  independent  simulations  was  performed  using  SARMAP  to  predict  the 
daylong  trajectories  at  successive  locations  along  the  paths  of  the  seven  SLDMB  in  the  Mid 
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Atlantic  Bight.  Comparison  of  the  model  predicted  and  observed  locations  at  the  end  of  one  day 
showed  a  median  separation  distance  of  7  km,  an  averaged  distance  of  9.3  km,  a  maximum 
distance  of  30  km  and  the  most  frequent  occurrence  distance  of  3  to  4  km.  The  error  in  location 
was  roughly  comparable  (about  80  percent)  to  the  distance  the  drifter  traveled  in  1  day  (5  to  18 
km/day).  Errors  in  predictions  were  significantly  higher  in  areas  where  the  CODAR  data  return 
rates  were  less  than  50  percent  (immediately  south  of  Long  Island  and  near  the  shelf  break). 
Estimates  of  the  dispersion  coefficients  necessary  as  input  to  the  model,  to  ensure  that  the 
predicted  and  observed  locations  were  within  the  model  predicted  search  area,  were  made  for 
each  simulation.  The  resulting  dispersion  coefficients  ranged  from  20  to  2,000  m2/sec,  with  most 
in  the  range  of  20  to  300  m2/sec.  The  median  value  was  about  90  m2/sec  and  the  most  frequent 
values  40  to  60  m2/sec.  The  largest  dispersion  coefficients  were  associated  with  areas  where  the 
CODAR  data  returns  were  lowest  and  hence,  not  representative. 

Estimates  of  the  median  horizontal  dispersion  coefficient  based  on  cluster  analysis,  CODAR 
velocity  errors,  and  SARMAP  trajectory  analysis  give  values  of  40  m2/sec,  30  to  60  m2/sec,  and 
(no  trajectory  analysis  in  Block  Island  Sound),  respectively  for  Block  Island  Sound.  Similar 
values  for  the  Mid  Atlantic  Bight  were  85,  200  to  400,  and  90  m2/sec.  This  range  of  values  is 
consistent  with  Okubo’s  4/3 rds  power  law,  which  gives  values  in  the  range  of  20  to  200  m2/sec, 
for  time  scales  of  5  to  30  days. 
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